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Abstract 

Harmonic inversion is introduced as a powerful tool for both the analysis of quan- 
tum spectra and semiclassical periodic orbit quantization. The method allows to 
circumvent the uncertainty principle of the conventional Fourier transform and to 
extract dynamical information from quantum spectra which has been unattainable 
before, such as bifurcations of orbits, the uncovering of hidden ghost orbits in com- 
plex phase space, and the direct observation of symmetry breaking effects. The 
method also solves the fundamental convergence problems in semiclassical periodic 
orbit theories - for both the Berry- Tabor formula and Gutzwiller's trace formula - 
and can therefore be applied as a novel technique for periodic orbit quantization, 
i.e., to calculate semiclassical eigenenergies from a finite set of classical periodic 
orbits. The advantage of periodic orbit quantization by harmonic inversion is the 
universality and wide applicability of the method, which will be demonstrated in 
this work for various open and bound systems with underlying regular, chaotic, 
and even mixed classical dynamics. The efficiency of the method is increased, i.e., 
the number of orbits required for periodic orbit quantization is reduced, when the 
harmonic inversion technique is generalized to the analysis of cross-correlated pe- 
riodic orbit sums. The method provides not only the eigenenergies and resonances 
of systems but also allows the semiclassical calculation of diagonal matrix elements 
and, e.g., for atoms in external fields, individual non-diagonal transition strengths. 
Furthermore, it is possible to include higher order terms of the ft expanded periodic 
orbit sum to obtain semiclassical spectra beyond the Gutzwiller and Berry- Tabor 
approximation . 

PACS numbers: 05.45. +b, 03.65.Sq 
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1 Introduction 



1.1 Motivation of semiclassical concepts 



Since the development of quantum mechanics in the early decades of this 
century quantum mechanical methods and computational techniques have be- 
come a powerful tool for accurate numerical calculations in atomic and molec- 
ular systems. The excellent agreement between experimental measurements 
and quantum calculations has silenced any serious critics on the fundamental 
concepts of quantum mechanics, and there are nowadays no doubts that quan- 
tum mechanics is the correct theory for microscopic systems. Nevertheless, 
there has been an increasing interest in semiclassical theories during recent 
years. The reasons for the resurgence of semiclassics are the following. Firstly, 
the quantum mechanical methods for solving multidimensional, non-integrable 
systems generically imply intense numerical calculations, e.g., the diagonaliza- 
tion of a large, but truncated Hamiltonian in a suitably chosen basis. Such 
calculations provide little insight into the underlying dynamics of the system. 
By contrast, semiclassical methods open the way to a deeper understanding of 
the system, and therefore can serve for the interpretation of experimental or 
numerically calculated quantum mechanical data in physical terms. Secondly, 
the relation between the quantum mechanics of microscopic systems and the 
classical mechanics of the macroscopic world is of fundamental interest and 
importance for a deeper understanding of nature. 

This relation was evident in the early days of quantum mechanics, when 
semiclassical techniques provided the only quantization rules, i.e., the WKB 
quantization of one-dimensional systems or the generalization to the Einstein- 
Brillouin-Keller (EBK) quantization [1-3] for systems with n degrees of free- 
dom. However, the EBK torus quantization is limited to integrable or at least 
near-integrable systems. In non-integrable systems the KAM tori are destroyed 
[4,5], a complete set of classical constants of motion does not exist any more, 
and therefore the eigenstates of the quantized system cannot be characterized 
by a complete set of quantum numbers. The "breakdown" of the semiclassical 
quantization rules for non- regular, i.e., chaotic systems was already recognized 
by Einstein in 1917 [1]. The failure of the "old" quantum mechanics to describe 
more complicated systems such as the helium atom [6], and, at the same time, 
the development and success of the "modern" wave mechanics are the reasons 
for little interest in semiclassical theories for several decades. The connection 
between wave mechanics and classical dynamics especially for chaotic systems 
remained an open question during that period. 



1.1.1 Basic semiclassical theories 

Chaotic systems: Gutzwiller's trace formula The problem was recon- 
sidered by Gutzwiller around 1970 [7,8]. Although a semiclassical quantiza- 
tion of individual eigenstates is, in principle, impossible for chaotic systems, 
Gutzwiller derived a semiclassical formula for the density of states as a whole. 
Starting from the exact quantum expression given as the trace of the Green's 
operator he replaced the exact Green's function with its semiclassical approx- 
imation. Applying stationary phase approximations he could finally write the 
semiclassical density of states as the sum of a smooth part, i.e., the Weyl 
term, and an oscillating sum over all periodic orbits of the corresponding clas- 
sical system. For this reason, Gutzwiller's theory is also commonly known 
as periodic orbit theory. Gutzwiller's semiclassical trace formula is valid for 
isolated and unstable periodic orbits, i.e., for fully chaotic systems with a 
complete hyperbolic dynamics. Examples of hyperbolic systems are n-disk re- 
pellers, as models for hyperbolic scattering [9-15], the stadium billiard [16-18], 
the anisotropic Kepler problem [8], and the hydrogen atom in magnetic fields 
at very high energies far above the ionization threshold [19-21]. Gutzwiller's 
trace formula is exact only in exceptional cases, e.g., for the geodesic flow on 
a surface with constant negative curvature [8]. In general, the semiclassical 
periodic orbit sum is just the leading order term of an infinite series in powers 
of the Planck constant, %. Methods to derive the higher order contributions 
of the h expansion are presented in Refs. [22-24]. 



Integrable systems: The Berry- Tabor formula For integrable systems 
a semiclassical trace formula was derived by Berry and Tabor [25,26]. The 
Berry- Tabor formula describes the density of states in terms of the periodic 
orbits of the system and is therefore the analogue of Gutzwiller's trace formula 
for integrable systems. The equation is known to be formally equivalent to the 
EBK torus quantization. A generalization of the Berry- Tabor formula to near- 
integrable systems is given in Refs. [27,28]. 



Mixed systems: Uniform semiclassical approximations Physical sys- 
tems are usually neither integrable nor exhibit complete hyperbolic dynamics. 
Generic systems are of mixed type, characterized by the coexistence of regu- 
lar torus structures and stochastic regions in the classical Poincare surface of 
section. A typical example is the hydrogen atom in a magnetic field [29-31], 
which undergoes a transition from near-integrable at low energies to complete 
hyperbolic dynamics at high excitation. In mixed systems the classical trajec- 
tories, i.e., the periodic orbits, undergo bifurcations. At bifurcations periodic 
orbits change their structure and stability properties, new orbits are born, or 
orbits vanish. A systematic classification of the various types of bifurcations is 
possible by application of normal form theory [32-35] , where the phase space 



structure around the bifurcation point is analyzed with the help of classical 
perturbation theory in local coordinates defined parallel and perpendicular to 
the periodic orbit. At bifurcation points periodic orbits are neither isolated 
nor belong to a regular torus in phase space. As a consequence the peri- 
odic orbit amplitudes diverge in both semiclassical expressions, Gutzwiller's 
trace formula and the Berry- Tabor formula, i.e., both formulae are not valid 
near bifurcations. The correct semiclassical solutions must simultaneously ac- 
count for all periodic orbits which participate at the bifurcation, including 
"ghost" orbits [36,37] in the complex generalization of phase space, which can 
be important near bifurcations. Such solutions are called uniform semiclassical 
approximations and can be constructed by application of catastrophe theory 
[38,39] in terms of catastrophe diffraction integrals. Uniform semiclassical ap- 
proximations have been derived in Refs. [40,41] for the simplest and generic 
types of bifurcations and in Refs. [37,42-44] for nongeneric bifurcations with 
higher codimension and corank. 

With Gutzwiller's trace formula for isolated periodic orbits, the Berry- Tabor 
formula for regular tori, and the uniform semiclassical approximations for or- 
bits near bifurcations we have, in principle, the basic equations for the semi- 
classical investigation of all systems with regular, chaotic, and mixed regular- 
chaotic dynamics. There are, however, fundamental problems in practical ap- 
plications of these equations for the calculation of semiclassical spectra, and 
the development of techniques to overcome these problems has been the ob- 
jective of intense research during recent years. 



1.1.2 Convergence problems of the semiclassical trace formulae 

The most serious problem of Gutzwiller's trace formula is that for chaotic 
systems the periodic orbit sum does not converge in the physically interesting 
energy region, i.e., on and below the real axis, where the eigenstates of bound 
systems and the resonances of open systems are located, respectively. For 
chaotic systems the sum of the absolute values of the periodic orbit terms 
diverges exponentially because of the exponential proliferation of the number 
of orbits with increasing periods. The convergence problems are similar for 
the quantization of regular systems with the Berry- Tabor formula, however, 
the divergence is algebraic instead of exponential in this case. 

Because of the technical problems encountered in trying to extract eigenvalues 
directly from the periodic orbit sum, Gutzwiller's trace formula has been used 
in many applications mainly for the interpretation of experimental or theoret- 
ical spectra of chaotic systems. For example the periodic orbit theory served 
as basis for the development of scaled- energy spectroscopy, a method where 
the classical dynamics of a system is fixed within long-range quantum spectra 
[45,46]. The Fourier transforms of the scaled spectra exhibit peaks at posi- 



tions given by the classical action of the periodic orbits. The action spectra 
can therefore be interpreted directly in terms of the periodic orbits of the sys- 
tem, and the technique of scaled-energy spectroscopy is now well established, 
e.g., for investigations of atoms in external fields [47-56]. The resolution of 
the Fourier transform action spectra is restricted by the uncertainty principle 
of the Fourier transform, i.e., the method allows the identification of usually 
short orbits in the low-dense part of the action spectra. A fully resolved ac- 
tion spectrum would require an infinite length of the original scaled quantum 
spectrum. 

Although periodic orbit theory has been very successful in the interpretation 
of quantum spectra, the extraction of individual eigenstates and resonances 
directly from periodic orbit quantization remains of fundamental interest and 
importance. As mentioned above the semiclassical trace formulae diverge at 
the physical interesting regions. They are convergent, however, at complex 
energies with positive imaginary part above a certain value, viz. the entropy 
barrier. Thus the problem of revealing the semiclassical eigenenergies and res- 
onances is closely related to finding an analytic continuation of the trace for- 
mulae to the region at and below the real axis. Several refinements have been 
introduced in recent years in order to transform the periodic orbit sum in the 
physical domain of interest to a conditionally convergent series by, e.g., us- 
ing symbolic dynamics and the cycle expansion [9,57,58], the Riemann-Siegel 
look-alike formula and pseudo-orbit expansions [59-61], surface of section tech- 
niques [62,63], and heat-kernel regularization [64-66]. These techniques are 
mostly designed for systems with special properties, e.g., the cycle expansion 
requires the existence and knowledge of a symbolic code and is most successful 
for open systems, while the Riemann-Siegel look-alike formula and the heat- 
kernel regularization are restricted to bound systems. Until now there is no 
universal method, which allows periodic orbit quantization of a large variety 
of bound and open systems with an underlying regular, mixed, or chaotic 
classical dynamics. 



1.2 Objective of this work 



The main objective of this work is the development of novel methods for (a) 
the analysis of quantum spectra and (6) periodic orbit quantization. The con- 
ventional action or recurrence spectra obtained by Fourier transformation of 
finite range quantum spectra suffer from the fundamental resolution problem 
because of the uncertainty principle of the Fourier transformation. The broad- 
ening of recurrence peaks usually prevents the detailed analysis of structures 
near bifurcations of classical orbits. We will present a method which allows, 
e.g., the detailed analysis of bifurcations and symmetry breaking and the study 
of higher order h corrections to the periodic orbit sum. For periodic orbit quan- 



tization the aim is the development of a universal method for periodic orbit 
quantization which does not depend on special properties of the system, such 
as the existence of a symbolic code or the applicability of a functional equa- 
tion. As will be shown both problems can be solved by applying methods for 
high resolution spectral analysis [67]. The computational techniques for high 
resolution spectral analysis have recently been significantly improved [68-70], 
and we will demonstrate that the state of the art methods are a powerful tool 
for the semiclassical quantization and analysis of dynamical systems. 



1.2.1 High precision analysis of quantum spectra 

Within the semiclassical approximation of Gutzwiller's trace formula or the 
Berry- Tabor formula the density of states of a scaled quantum spectrum is 
a superposition of sinusoidal modulations as a function of the energy or an 
appropriate scaling parameter. The frequencies of the modulations are deter- 
mined by the periods, i.e., the classical action of orbits. The amplitudes and 
phases of the oscillations depend on the stability properties of the periodic 
orbits and the Maslov indices. To extract information about the classical dy- 
namics from the quantum spectrum it is therefore quite natural to Fourier 
transform the spectrum from "energy" domain to "time" (or, for scaled spec- 
tra, "action") domain [45-47,71-77]. The periodic orbits appear as peaks at 
positions given by the periods of orbits, and the peak heights exhibit the ampli- 
tudes of periodic orbit contributions in the semiclassical trace formulae. The 
comparison with classical calculations allows the interpretation of quantum 
spectra in terms of the periodic orbits of the corresponding classical system. 
However, the analyzed experimental or theoretical quantum spectra are usu- 
ally of finite length, and thus the sinusoidal modulations are truncated when 
Fourier transformed. The truncation implies a broadening of peaks in the time 
domain spectra because of the uncertainty principle of the Fourier transform. 
The widths of the recurrence peaks are determined by the total length of the 
original spectrum. The uncertainty principle prevents until now precision tests 
of the semiclassical theories because neither the peak positions (periods of or- 
bits) nor the amplitudes can be obtained from the time domain spectra to 
high accuracy. Furthermore, the uncertainty principle implies an overlapping 
of broadened recurrence peaks when the separation of neighboring periods is 
less than the peak widths. In this case individual periodic orbit contributions 
cannot be revealed any more. This is especially a disadvantage, e.g., when 
following, in quantum spectra of systems with mixed regular-chaotic dynam- 
ics, the bifurcation tree of periodic orbits. All orbits involved in a bifurcation, 
including complex "ghost" orbits [36,37,44], have nearly the same period close 
to the bifurcation point. The details of the bifurcations, especially of catastro- 
phes with higher codimension or corank [38,39], cannot be resolved with the 
help of the conventional Fourier transform. The same is true for the effects of 
symmetry breaking, e.g., breaking of the cylindrical symmetry of the hydrogen 



atom in crossed magnetic and electric fields [54] or the "temporal symmetry 
breaking" [55,56] of atoms in oscillating fields. The symmetry breaking should 
result in a splitting of peaks in the Fourier transform recurrence spectra. Until 
now this phenomenon could only be observed indirectly via constructive and 
destructive interference effects of orbits with broken symmetry [54,55]. 

To overcome the resolution problem in the conventional recurrence spectra and 
to achieve a high resolution analysis of quantum spectra it is first of all neces- 
sary to mention that the uncertainty principle of the Fourier transformation is 
not a fundamental restriction to the resolution of the time domain recurrence 
spectra, and is not comparable to the fundamental Heisenberg uncertainty 
principle in quantum mechanics. This can be illustrated for the simple ex- 
ample of a sinusoidal function with only one frequency. The knowledge of 
the signal at only three points is sufficient to recover, although not uniquely, 
the three unknown parameters, viz. the frequency, amplitude and phase of 
the modulation. Similarly, a superposition of N sinusoidal functions can be 
recovered, in principle, from a set of at least 3N data points (£«,/«)• If the 
data points are exact the spectral parameters {u>k, A k , <pk} can also be exactly 
determined from a set of 3N nonlinear equations 

N 

fi = ^A k sm(u k ti-(t>k) , (1.1) 

fc=i 

i.e., there is no uncertainty principle. The main difference between this proce- 
dure and the Fourier transformation is that we use the linear superposition of 
sinusoidal functions as an ansatz for the functional form of our data, and the 
frequencies, amplitudes, and phases are adjusted to obtain the best approxi- 
mation to the data points. By contrast, the Fourier transform is not based on 
any constraints. In case of a discrete Fourier transform (DFT) the frequencies 
are chosen on a (usually equidistant) grid, and the amplitudes are determined 
from a linear set of equations, which can be solved numerically very efficiently, 
e.g., by fast Fourier transform (FFT). 

As mentioned above the high precision spectral analysis requires the numerical 
solution of a nonlinear set of equations, and the availability of efficient algo- 
rithms has been the bottleneck for practical applications in the past. Several 
techniques have been developed [67], however, most of them are restricted - 
for reasons of storage requirements, computational time, and stability of the 
algorithm - to signals with quite low number of frequencies. By contrast, the 
number of frequencies (periodic orbits) in quantum spectra of chaotic systems 
is infinite or at least very large, and applied methods for high precision spec- 
tral analysis must be able to handle signals with large numbers of frequencies. 
This requirement is fulfilled by the method of harmonic inversion by filter- 
diagonalization, which was recently developed by Wall and Neuhauser [68] and 
significantly improved by Mandelshtam and Taylor [69,70]. The decisive step 



is to recast the nonlinear set of equations as a generalized eigenvalue problem. 
Using an appropriate basis set the generalized eigenvalue equation is solved 
with the filter-diagonalization method, which means that the frequencies of 
the signal can be determined within small frequency windows and only small 
matrices must be diagonalized numerically even when the signal is composed 
of a very large number of sinusoidal oscillations. 

We will introduce harmonic inversion as a powerful tool for the high precision 
analysis of quantum spectra, which allows to circumvent the uncertainty prin- 
ciple of the conventional Fourier transform analysis and to extract previously 
unattainable information from the spectra. In particular the following items 
are investigated: 

• High precision check of semiclassical theories: The analysis of spectra allows 
a direct quantitative comparison of the quantum and classical periodic orbit 
quantities to many significant digits. 

• Uncovering of periodic orbit bifurcations in quantum spectra, and the veri- 
fication of ghost orbits and uniform semiclassical approximations. 

• Direct observation of symmetry breaking effects. 

• Quantitative interpretation of the differences between quantum and semi- 
classical spectra in terms of the % expansion of the periodic orbit sum. 

Results will be presented for various physical systems, e.g., the hydrogen atom 
in external fields, the circle billiard, and the three disk scattering problem. 



1.2.2 Periodic orbit quantization 

The analysis of quantum spectra by harmonic inversion and, e.g., the observa- 
tion of "ghost" orbits, symmetry breaking effects, or higher order H corrections 
to the periodic orbit contributions provides a deeper understanding of the re- 
lation between quantum mechanics and the underlying classical dynamics of 
the system. However, the inverse procedure, i.e., the calculation of semiclas- 
sical eigenenergies and resonances directly from classical input data is of at 
least the same importance and even more challenging. As mentioned above 
the periodic orbit sums suffer from fundamental convergence problems in the 
physically interesting domain and much effort has been undertaken in recent 
years to overcome these problems [9,57-66]. Although many of the refinements 
which have been introduced are very efficient for a specific model system or a 
class of systems, they all suffer from the disadvantage of non-universality The 
cycle expansion technique [9,57,58] requires a completely hyperbolic dynam- 
ics and the existence of a symbolic code. The method is most efficient only 
for open systems, e.g., for three-disk or n-disk pinball scattering [9-15]. By 
contrast, the Riemann-Siegel look-alike formula and pseudo orbit expansion of 
Berry and Keating [59-61] can only be applied for bound systems. The same 



is true for surface of section techniques [62,63] and heat-kernel regularization 
[64-66]. We will introduce high resolution spectral analysis, and in particular 
harmonic inversion by filter-diagonalization as a novel and universal method 
for periodic orbit quantization. Universality means that the method can be ap- 
plied to open and bound systems with regular and chaotic classical dynamics 
as well. 

Formally the semiclassical density of states (more precisely the semiclassical 
response function) can be written as the Fourier transform of the periodic 
orbit recurrence signal 

C sc (s) = J2Ao5(s-s po ), (1.2) 

po 



with A po and s po the periodic orbit amplitudes and periods (actions), respec- 
tively. If all orbits up to infinite length are considered the Fourier transform 
of C sc (s) is equivalent to the non-convergent periodic orbit sum. For chaotic 
systems a numerical search for all periodic orbits is impossible and, further- 
more, does not make sense because the periodic orbit sum does not converge 
anyway. On the other hand, the truncation of the Fourier integral at finite 
maximum period s max yields a smoothed spectrum only [78]. However, the 
low resolution property of the spectrum can be interpreted as a consequence 
of the uncertainty principle of the Fourier transform. We can now argue in an 
analogous way as in the previous section that the uncertainty principle can be 
circumvented with the help of high resolution spectral analysis, and propose 
the following procedure for periodic orbit quantization by harmonic inversion. 

Let us assume that the periodic orbits with periods < s < s max are available 
and the semiclassical recurrence function C sc (s) has been constructed. This 
signal can now be harmonically inverted and thus adjusted to the functional 
form of the quantum recurrence function C qm (s), which is a superposition 
of sinusoidal oscillations with frequencies given by the quantum mechanical 
eigenvalues of the system. The frequencies obtained by harmonic inversion 
of C sc (s) should therefore be the semiclassical approximations to the exact 
eigenenergies. The universality and wide applicability of this novel quantiza- 
tion scheme follows from the fact that the periodic orbit recurrence signal 
C sc (s) can be obtained for a large variety of systems because the calcula- 
tion does not depend on special properties of the system, such as boundness, 
ergodicity, or the existence of a symbolic dynamics. For systems with underly- 
ing regular or chaotic classical dynamics the amplitudes A po of periodic orbit 
contributions are directly obtained from the Berry- Tabor formula [25,26] and 
Gutzwiller's trace formula [8], respectively. 

As mentioned above the harmonic inversion technique requires the knowledge 
of the signal C sc (s) up to a finite maximum period s max . The efficiency of the 
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quantization method strongly depends on the signal length which is required 
to obtain a certain number of eigenenergies. In chaotic systems periodic orbits 
proliferate exponentially with increasing period and usually the orbits must 
be searched numerically. It is therefore highly desirable to use the shortest 
possible signal for periodic orbit quantization. The efficiency of the method 
can be improved if not just the single signal C sc (s) is harmonically inverted 
but additional classical information obtained from a set of smooth and linearly 
independent observables is used to construct a semiclassical cross-correlated 
periodic orbit signal. The cross-correlation function can be analyzed with a 
generalized harmonic inversion technique and allows calculating semiclassi- 
cal eigenenergies from a significantly reduced set of orbits or alternatively to 
improve the accuracy of spectra obtained with the same set of orbits. 



However, the semiclassical eigenenergies deviate - apart from a few exceptions, 
e.g., the geodesic flow on a surface with constant negative curvature [8] - from 
the exact quantum mechanical eigenvalues. The reason is that Gutzwiller's 
trace formula and the Berry- Tabor formula are only the leading order terms 
of the h expansion of periodic orbit contributions [22-24]. It will be shown how 
the higher order h corrections of the periodic orbit sum can be used within 
the harmonic inversion procedure to improve, order by order, the semiclassical 
accuracy of eigenenergies, i.e., to obtain eigenvalues beyond the Gutzwiller and 
Berry- Tabor approximation. 



1.3 Outline 



The manuscript is organized as follows. In Chapter 2 the high precision anal- 
ysis of quantum spectra is discussed. After general remarks on Fourier trans- 
form recurrence spectra in Section 2.1 we introduce in Section 2.2 harmonic 
inversion as a tool to circumvent the uncertainty principle of the conventional 
Fourier transformation [79]. In Section 2.3 a precision check of the periodic 
orbit theory is demonstrated by way of example of the hydrogen atom in a 
magnetic field. For the hydrogen atom in external fields we furthermore illus- 
trate that dynamical information which has been unattainable before can be 
extracted from the quantum spectra. In Section 2.4 we investigate in detail 
the quantum manifestations of bifurcations of orbits related to a hyperbolic 
umbilic catastrophe [44,80] and a butterfly catastrophe [37,81], and in Section 
2.5 we directly uncover effects of symmetry breaking in the quantum spectra. 
In Section 2.6 we analyze by way of example of the circle billiard and the three 
disk system the difference between the exact quantum and the semiclassical 
spectra. It is shown that the deviations between the spectra can be quantita- 
tively interpreted in terms of the next order h contributions of the % expanded 
periodic orbit sum. 
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In Chapter 3 we propose harmonic inversion as a novel technique for peri- 
odic orbit quantization [82,83]. The method is introduced in Section 3.1 for a 
mathematical model, viz. the calculation of zeros of Riemann's zeta function. 
This model shows a formal analogy with the semiclassical trace formulae, and 
is chosen for the reasons that no extensive numerical periodic orbit search is 
necessary and the results can be directly compared to the exact zeros of the 
Riemann zeta function. In Section 3.2 the method is derived for the periodic 
orbit quantization of physical systems and is demonstrated in Section 3.3 for 
the example of the three disk scattering system with fully chaotic (hyperbolic) 
classical dynamics [82,83], and in Section 3.4 for the hydrogen atom in a mag- 
netic field as a prototype example of a system with mixed regular-chaotic 
dynamics [84]. In Section 3.5 the efficiency of the method is improved by a 
generalization of the technique to the harmonic inversion of cross-correlated 
periodic orbit sums, which allows to significantly reduce the number of or- 
bits required for the semiclassical quantization [85], and in Section 3.6 we 
derive the concept for periodic orbit quantization beyond the Gutzwiller and 
Berry- Tabor approximation by harmonic inversion of the H expansion of the 
periodic orbit sum [86]. The methods of Sections 3.5 and 3.6 are illustrated in 
Section 3.7 for the example of the circle billiard [85-88]. Finally, in Section 3.8 
we demonstrate the semiclassical calculation of individual transition matrix 
elements for atoms in external fields [89]. 

Chapter 4 concludes with a summary and outlines possible future applications, 
e.g., the analysis of experimental spectra and the periodic orbit quantization 
of systems without scaling properties. 

Computational details of harmonic inversion by filter-diagonalization are pre- 
sented in Appendix A, and the calculation and asymptotic behavior of some 
catastrophe diffraction integrals is discussed in Appendix C. 



2 High precision analysis of quantum spectra 



Semiclassical periodic orbit theory [7,8] and closed orbit theory [90-93] have 
become the key for the interpretation of quantum spectra of classically chaotic 
systems. The semiclassical spectra at least in low resolution are given as the 
sum of a smooth background and a superposition of modulations whose ampli- 
tudes, frequencies, and phases are solely determined by the closed or periodic 
orbits of the classical system. For the interpretation of quantum spectra in 
terms of classical orbits it is therefore most natural to obtain the recurrence 
spectrum by Fourier transforming the energy spectrum to the time domain. 
Each closed or periodic orbit should show up as a sharp <5-peak at the re- 
currence time (period), provided, first, the classical recurrence times do not 
change within the whole range of the spectrum and, second, the Fourier inte- 
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gral is calculated along an infinite energy range. Both conditions are usually 
not fulfilled. However, the problem regarding the energy dependence of recur- 
rence times can be solved in systems possessing a classical scaling property 
by application of scaling techniques. The second condition is never fulfilled 
in practice, i.e., the length of quantum spectra is always restricted either by 
experimental limitations or, in theoretical calculations, by the growing dimen- 
sion of the Hamiltonian matrix which has to be diagonalized numerically. The 
given length of a quantum spectrum determines the resolution of the quantum 
recurrence spectrum due to the uncertainty principle, AE ■ AT ~ h, when the 
conventional Fourier transform is used. Only those closed or periodic orbits 
can be clearly identified quantum mechanically which appear as isolated non- 
overlapping peaks in the quantum recurrence spectra. This is especially not 
the case for orbits which undergo bifurcations at energies close to the bifurca- 
tion point. As will be shown the resolution of quantum recurrence spectra can 
be significantly improved beyond the limitations of the uncertainty principle 
of the Fourier transform by methods of high resolution spectral analysis. 

In the following we first review the conventional analysis of quantum spectra 
by Fourier transformation and discuss on the example of the hydrogen atom 
in a magnetic field the achievements and limitations of the Fourier transform 
recurrence spectra. In Section 2.2 we review methods of high resolution spec- 
tral analysis which can serve to circumvent the uncertainty principle of the 
Fourier transform. In Sections 2.3 to 2.6 the harmonic inversion technique is 
applied to calculate high resolution recurrence spectra beyond the limitations 
of the uncertainty principle of the Fourier transformation from experimental 
or theoretical quantum spectra of finite length. The method allows to reveal in- 
formation about the dynamics of the system which is completely hidden in the 
Fourier transform recurrence spectra. In particular, it allows to identify real 
orbits with nearly degenerate periods, to detect complex "ghost" orbits which 
are of importance in the vicinity of bifurcations [36,37,44], and to investigate 
higher order % corrections of the periodic orbit contributions [22-24]. 



2. 1 Fourier transform recurrence spectra 



According to periodic orbit theory [7,8] the semiclassical density of states can 
be written as the sum of a smooth background Qo(E) and oscillatory modula- 
tions induced by the periodic orbits, 

-1 oo 

g(E) = g (E) + -Re £ £ A po ,{E)e^^' h , (2.1) 

" po r=l 



with A p0yr (E) the amplitudes of the modulations and the classical actions 
Sp (E) of a primitive periodic orbit (po) determining the frequencies of the 
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oscillations. Linearizing the action around E = Eq yields 

dS 



S(E) « S(E ) + (]e 



(E - E ) = S(E ) + T (E - E ) , (2.2) 

E 



with To the time period of the orbit at energy Eq. When Eq. (2.2) is inserted 
in (2.1) the semiclassical density of states is locally given as a superposition of 
sinusoidal oscillations, and it might appear that the problem of identifying the 
amplitudes A po and time periods T po that contribute to the quantum spec- 
trum g(E) can be solved by Fourier transforming g(E) to the time domain, i.e., 
each periodic orbit should show up as a 5-peak at the recurrence time t — T . 
However, fully resolved recurrence spectra can be obtained only if the two fol- 
lowing conditions are fulfilled. First, the Fourier integral is calculated along an 
infinite energy range, and, second, the classical recurrence times do not change 
within the whole range of the spectrum. The first condition is usually not ful- 
filled when quantum spectra are obtained from an experimental measurement 
or a numerical quantum calculation. In that case the Fourier transformation 
is restricted to a finite energy range and the resolution of the time domain 
spectrum is limited by the uncertainty principle of the Fourier transform. It 
is the main objective of this Chapter to introduce high resolution methods for 
the spectral analysis which allow to circumvent the uncertainty principle of 
the Fourier transform and to obtain fully resolved recurrence spectra from the 
analysis of quantum spectra with finite length. However, the second condition 
is usually not fulfilled either. Both the amplitudes and recurrence times of 
periodic orbits are in general nontrivial functions of the energy E, and, even 
worse, the whole phase space structure changes with energy when periodic or- 
bits undergo bifurcations. For the interpretation of quantum spectra in terms 
of periodic orbit quantities it is therefore in general not appropriate to analyze 
the frequencies of long ranged spectra g(E). The problems due to the energy 
dependence of periodic orbit quantities have been solved by the development 
and application of scaling techniques. 



Scaling techniques Many systems possess a classical scaling property in 
that the classical dynamics does not depend on an external scaling parameter 
w and the action of trajectories varies linearly with w. Examples are systems 
with homogeneous potentials, billiard systems, or the hydrogen atom in exter- 
nal fields. In billiard systems the shapes of periodic orbits are solely determined 
by the geometry of the borders, and the classical action depends on the length 
L po of the periodic orbit, S po = hkL po , with k the wave number. For a particle 
with mass m moving in a billiard system it is therefore most appropriate to 
take the wave number as the scaling parameter, i.e., w = k = y/2mE/h. For a 
system with a homogeneous potential V(q) with V(cq) = c a V(q) only the size 
of periodic orbits but not their shape changes with varying energy E. Introduc- 
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ing a scaling parameter as a power of the energy, w = E 1 / 01 '*" 1 / 2 , the classical 
action of a periodic orbit is obtained as S po (E) = (E / E o y/ a+1 / 2 S po (E ) , with 
Eq being a reference energy, i.e., the action depends linearly on a scaling pa- 
rameter defined as w = (E/Eq) 1 /^ 1 / 2 . For example the Coulomb potential is 
a homogeneous potential with a = — 1 and the bound Coulomb spectrum at 
negative energies (Eq = —1/2) is transformed by the scaling procedure to a 
simple spectrum with equidistant lines at w = n with n = 1,2,... the princi- 
pal quantum number. For atoms in external magnetic and electric fields the 
shape of periodic orbits changes if the field strengths are fixed and only the 
energy is varied. However, these systems possess a scaling property if both 
the energy and the field strengths are simultaneously varied. Details for the 
hydrogen atom in a magnetic field will be given below. The scaling parameter 
plays the role of an inverse effective Planck constant, w = h~g. In theoretical 
investigations it is even possible to apply scaling techniques to general sys- 
tems with non-homogeneous potentials if the Planck constant is formally used 
as a variable parameter. Non-homogeneous potentials are important, e.g., in 
molecular dynamics, and the genaralized scaling technique has been applied 
in Ref. [94] to analyze quantum spectra of the HO2 molecule. 

When the scaling technique is applied to a quantum system the scaling param- 
eter w is quantized, i.e., bound systems exhibit sharp lines at real eigenvalues 
w = Wk and open systems have resonances related to complex poles w = Wk 
of the scaled Green's function G + (w). By varying the scaled energy a direct 
comparison of the quantum recurrence spectra with the bifurcation diagram 
of the underlying classical system is possible [45,47]. 

The semiclassical approximation to the scaled spectrum is given by 

-1 00 

g(w) = qo(w) + -Re £ £ A p0tl e irs »° w , (2.3) 

™ po r=l 



with 



■5po 



j> prfq = — f pdq (2.4) 



po po 



the scaled action of a primitive periodic orbit. The vectors q and p are the 
coordinates and momenta of the scaled Hamiltonian. In contrast to Eq. 2.1 
the periodic orbit sum (2.3) is a superposition of sinusoidal oscillations as a 
function of the scaling parameter w. Therefore the scaled spectra g(w) can 
be Fourier transformed along arbitrarily long ranges of w to generate Fourier 
transform recurrence spectra of in principle arbitrarily high resolution, i.e., 
yielding sharp 5-peaks at the positions of the scaled action s = s po of peri- 
odic orbits. The high resolution analysis of quantum spectra in the following 
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sections is possible only in conjunction with the application of the scaling 
technique. 

The periodic orbit amplitudes A por in Gutzwiller's trace formula (2.3) for 
scaled systems are given by 



A po , = / po er-W , (2.5) 

v /|det(M; o -/)| 

with M po and /x po the monodromy matrix and Maslov index of the primi- 
tive periodic orbits, respectively, and r the repetition number of orbits. The 
monodromy matrix M is the stability matrix restricted to deviations perpen- 
dicular to a periodic orbit after period time T. We here discuss systems with 
two degrees of freedom. If Sq(0) is a small deviation perpendicular to the orbit 
in coordinate space at time t — and Sp(0) an initial deviation in momen- 
tum space, the corresponding deviations at time t = T are related to the 
monodromy matrix [93,95]: 

f^»)=Mf^)) = f'"» m -)(^)). (2.6) 

\5p(T) J \op{0) J \m 2 i m 2 2 J \op{0) J 



To compute M one considers an initial deviation solely in coordinate space 
to obtain the matrix elements ran and TO21, and an initial deviation solely 
in momentum space to obtain 77112 and m22- In practice a linearized system 
of differential equations obtained by differentiating Hamilton's equations of 
motion with respect to the phase space coordinates is numerically integrated. 
The Maslov index /i po increases by one every time the trajectory passes a con- 
jugate point or a caustic. Therefore the amplitudes A po ^ T are complex numbers 
containing phase information determined by the Maslov indices of orbits. The 
classical actions s po are usually real numbers, although they can be complex 
in general. Non-real actions s po indicate "ghost" orbits [36,37,44] which exist 
in the complex continuation of the classical phase space. 

As mentioned above the problem of identifying A po and s po that contribute to 
the quantum spectrum g(w) can in principle be solved by Fourier transforming 
g(w) to the action domain, 

C(s) = — J [g(w) - g (w)} e~* sw dw . (2.7) 



If the quantum spectrum g(w) has infinite length and the periodic orbits (i.e., 
the actions s po ) are real the Fourier transform of Eq. 2.3 indeed results in a 
fully resolved recurrence spectrum 
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C(s) = ±- f [g(w)-g (w)]e-' tsw dw 

— oo 
1 oo 

= 7T- E E 0W(« - rs po ) + A^ r <J(s + rs po )) . (2.8) 

Z7r po r=l 

The periodic orbits are identified as <5-peaks in the recurrence spectrum at 
positions s = ±s po (and their repetitions at s = ±rs po ) and the recurrence 
strengths are given as the amplitudes A por (Eq. 2.5) of Gutzwiller's trace 
formula. However, for finite range spectra the 5-functions in Eq. 2.8 must be 
replaced apart from a phase factor with 

sin((u? 2 -iui)s/2) 
6{s) — ► . 

7TS 



Unfortunately, the recurrence peaks are broadened by the uncertainty principle 
of the Fourier transform and furthermore the function sin((u>2 — Wi)s/2)/tts 
has side-peaks which are not related to recurrences of periodic orbits but are 
solely an undesirable effect of the sharp cut of the Fourier integral at w± and W2- 
In complicated recurrence spectra it can be difficult or impossible to separate 
the physical recurrences and the unphysical side-peaks. The occurrence of 
side-peaks can be avoided by multiplying the spectra with a window function 
h(w), i.e. 

w 2 

C(s) = —J h(w) [g(w) - qo(w)] e~ lsw dw , (2.9) 



where h(w) is equal to one at the center of the spectrum and decreases 
smoothly close to zero at the borders u>i and W2 of the Fourier integral. An 
example is a Gaussian window 

h( \ I ( w - w o) 2 \ 



centered at wq = (wi + w-^jl and with sufficiently chosen width a ~ [wi — 
Wi)/6. The Gaussian window suppresses unphysical side-peaks in the recur- 
rence spectra, however, the decrease of the side-peaks is paid by an additional 
broadening of the central peak. Various other types of window functions h(w) 
have been used as a compromise between the least uncertainty of the central 
recurrence peak and the optimal suppression of side-peaks. However, the fun- 
damental problem of the uncertainty principle of the Fourier transform cannot 
be solved with any window function h(w). 
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As a first example for the analysis of quantum spectra we introduce the hydro- 
gen atom in a magnetic field (for reviews see [29-31]) given by the Hamiltonian 
[in atomic units, magnetic field strength 7 = 5/(2.35 x 10 5 T), angular mo- 
mentum L z = 0] 

tf=ip 2 -i + i 7 V + 2/ 2 )- (2.10) 

2 r 8 



As mentioned above this system possesses a scaling property. Introducing 
scaled coordinates, f = 7 2/,3 r, and momenta, p = 7 _1//3 p, the classical dy- 
namics of the scaled Hamiltonian 

H= l -tf-\+ l -^ + ?) = E 1 - 21 *. (2.11) 



does not depend on two parameters, the energy E and the magnetic field 
strength 7, but solely on the scaled energy 

E = E-f- 2/3 . (2.12) 

The classical action of the trajectories scales as 

S po = s P o7~ 1/3 = s po w . (2.13) 

Based on the scaling relations of the classical Hamiltonian the experimental 
technique of scaled energy spectroscopy was developed [45-47] . Experimental 
spectra on the hydrogen atom at constant scaled energy have been measured 
by varying the magnetic field strength linearly on a scale w = 7 -1 ' 3 , adjusting 
simultaneously the energy E (via the wave length of the exciting laser light) 
so that the scaled energy E = E'j -2 ' 3 is kept constant at a given value. 
The spectra have been Fourier transformed in the experimentally accessible 
range 36 < 7 -1 / 3 < 50. Experimental recurrence spectra at scaled energies 
—0.3 < E < are presented in Fig. la. The overlay exhibits a well-structured 
system of clustered branches of resonances in the scaled energy-action plane. 
For comparison Fig. lb presents the energy-action spectrum of the closed 
classical orbits, i.e. the scaled action of orbits as a function of the energy. As 
can be seen the clustered branches of resonances in the experimental recurrence 
spectra (Fig. la) well resemble the classical bifurcation tree of closed orbits in 
Fig. lb, although a one to one comparison between recurrence peaks and closed 
orbits is limited by the finite resolution of the quantum recurrence spectra. 
Closed orbits bifurcating from the /ith repetition of the orbit parallel to the 
magnetic field axis are called "vibrators" VT, and orbits bifurcating from the 
/ith repetition of the perpendicular orbit are called "rotators" R v in Fig. 1. 
Orbits X M are created mainly in tangent bifurcations. The graphs of some 



closed orbits are given in Fig. 2. A detailed comparison of the peak heights 
of resonances in the experimental recurrence spectra with the semiclassical 
amplitudes obtained from closed orbit theory can be found in Ref. [47]. 

Scaled energy spectroscopy has become a well established method to investi- 
gate the dynamics of hydrogenic and nonhydrogenic atoms in external fields. 
Recurrence spectra of helium in a magnetic field are presented in Ref. [48,49]. 
Recurrence peaks of the helium atom which cannot be identified by hydrogenic 
closed orbits have been explained by scattering of the highly excited electron 
at the ionic core [49,75,76]. Rubidium has been studied in crossed electric and 
magnetic fields in Ref. [50] and unidentified recurrence peaks have been inter- 
preted in terms of classical core scattering [77]. The Stark effect on lithium has 
been investigated by the MIT group [51-53]. For atoms in an electric field F 
the scaling parameter is w = F -1 / 4 , and EF~~ l l 2 is the scaled energy. Experi- 
mental recurrence spectra of lithium in an electric field and the corresponding 
closed orbits are presented in Fig. 3. Strong recurrence peaks occur close to 
the bifurcations of the parallel orbit and its repetitions marked by open circles 
in Fig. 3. The new orbits created in bifurcations have almost the same action 
as the corresponding return of the parallel orbit, and, similar as in Fig. 1 for 
the hydrogen atom in a magnetic field, the small splittings of recurrence peaks 
are not resolved in the experimental recurrence spectra. 

It is important to note that the finite resolution of the experimental recurrence 
spectra in Figs. 1 and 3 is not caused by the finite bandwidth of the exciting 
laser but results as discussed above from the finite length of the Fourier trans- 
formed scaled spectra and the uncertainty principle of the Fourier transform. 
This can be illustrated by analyzing a numerically computed spectrum instead 
of an experimentally measured spectrum. We study the hydrogen atom in a 
magnetic field at constant scaled energy E = —0.1. At this energy the classical 
dynamics is completely chaotic and all periodic orbits are unstable. We cal- 
culated 9715 states in the region w < 140 by numerical diagonalization of the 
Hamiltonian matrix in a complete basis set. For details of the quantum calcu- 
lations see, e.g., Ref. [96]. In Fig. 4 the quantum density of states is analyzed 
by the conventional Fourier transform. To get rid of unphysical side-peaks the 
spectrum was multiplied with a Gaussian function h(w) with width a chosen 
in accordance with the total length of the quantum spectrum. Fig. 4 clearly 
exhibits recurrence peaks which can be related to classical periodic orbits. 
However, the widths of the peaks are approximately As/2tt = 0.03, and it is 
impossible to determine the periods of the classical orbits to higher accuracy 
than about 0.01 from the Fourier analysis of the quantum spectrum. Recur- 
rence peaks of orbits with similar periods overlap, as can be clearly seen around 
s/2n = 2.1, and at least guessed, e.g., at s/2n = 1.1, or s/2n = 2.6. A precise 
determination of the amplitudes is impossible especially for the overlapping 
peaks. Furthermore, the Fourier transform does not allow to distinguish be- 
tween real and ghost orbits. In the following we will demonstrate that the 
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quality of recurrence spectra can be significantly improved by application of 
state of the art methods for high resolution spectral analysis. 



2.2 Circumventing the uncertainty principle 



Instead of using the standard Fourier analysis, to extract the amplitudes and 
actions we propose to apply methods for high resolution spectral analysis. 
The fundamental difference between the Fourier transform and high resolu- 
tion methods is the following. Assume a complex signal C(t) given on an 
equidistant grid of points t n = t + nAt as a superposition of exponential 
functions, i.e. 

C(t n ) = c n = J2 d k e~ UnUJk , n = 1, 2, . . . , N . (2.14) 

k 



In the Discrete Fourier Transform (DFT) N real frequencies 

2irk 
^ = iVAt' k = 1 > 2 >---> N ( 2 - 15 ) 



are fixed and evenly spaced. The N complex amplitudes dk of the Fourier 
transform are determined by solving a linear set of equations, yielding 

1 N 
d k = - £ c n e u ^ , k = l,2,...,N. (2.16) 



The sums in Eq. 2.16 can be calculated numerically very efficiently, e.g., by the 
Fast Fourier Transform algorithm (FFT). The resolution Auj of the Fourier 
transform is controlled by the total length NAt of the signal, 

2tt 
NAt y ' 



which is the spacing between the grid points in the frequency domain. Only 
those spectral features that are separated from each other by more than Aw 
can be resolved. This is referred to as the "uncertainty principle" of the Fourier 
transform. By contrast, both the amplitudes d\. and the frequencies u>k in Eq. 
2.14 are free parameters when methods for high resolution spectral analysis 
are applied. Because the frequencies u>k are free adjusting parameters they 
are allowed to appear very close to each other and therefore the resolution is 
practically infinite. Using the data points of the signal C(t) at t — t n the set of 
in general complex parameters {dk, ojk] are given as the solution of a nonlinear 
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set of equations. Unfortunately, the nonlinear set of equations does not have a 
solution in closed form similar to Eq. 2.16 for the discrete Fourier transform. In 
fact, the numerical calculation of the parameters {dk,uJk} from the nonlinear 
set of equations is the central and nontrivial problem of all methods for high 
resolution spectral analysis. 

The numerical harmonic inversion of a given signal like Eq. 2.14 is a fundamen- 
tal problem in physics, electrical engineering and many other diverse fields. 
The problem has already been addressed (in a slightly different form) in the 
18th century by Baron de Prony, who converted the nonlinear set of equations 
(2.14) to a linear algebra problem. There are several approaches related to the 
Prony method used for a high resolution spectral analysis of short time signals, 
such as the modern versions of the Prony method, MUSIC (Multiple Signal 
Classification) and ESPRIT (Estimation of Signal Parameters via Rotational 
Invariance Technique) [67,97,98]. As opposed to the Fourier transform these 
methods are not related to a linear (unitary) transformation of the signal and 
are highly nonlinear by their nature. However, the common feature present 
in most versions of these methods is converting the nonlinear fitting problem 
to a linear algebraic one. Note also that while ESPRIT uses exclusively lin- 
ear algebra, the Prony method or MUSIC require some additional search in 
the frequency domain which makes them less efficient. To our best knowledge 
none of these nonlinear methods is able to handle a signal (2.14) that con- 
tains "too many" frequencies as they lead to unfeasibly large and typically 
ill-conditioned linear algebra problems [70]. This is especially the case for the 
analysis of quantum spectra because such spectra cannot be treated as "short 
signals" and contain a high number of frequencies given by (see below) the 
number of periodic orbits of the underlying classical system. Decisive progress 
in the numerical techniques for harmonic inversion has recently been achieved 
by Wall and Neuhauser [68]. Their method is conceptually based on the orig- 
inal filter-diagonalization method of Neuhauser [99] designed to obtain the 
eigenspectrum of a Hamiltonian operator in any selected small energy range. 
The main idea is to associate the signal C(t) (Eq. 2.14) with an autocorrelation 
function of a suitable dynamical system, 

C(t) = f^o,e-^ ) , (2.18) 



where the brackets define a complex symmetric inner product (i.e., no com- 
plex conjugation). Eq. 2.18 establishes an equivalence between the problem 
of extracting information from the signal C(t) with the one of diagonalizing 
the evolution operator exp(— iCl) of the underlying dynamical system. The 
frequencies u>k of the signal C(t) are the eigenvalues of the operator Cl, i.e. 

ft|$ fe )=u; fe |$ fe ), (2.19) 
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and the amplitudes dj~ are formally given as 

4 = ($ ,$ fc ) 2 • (2.20) 



After introducing an appropriate basis set, the operator exp(— iCl) can be diag- 
onalized using the method of filter- diagonalization [68-70]. Operationally this 
is done by solving a small generalized eigenvalue problem whose eigenvalues 
yield the frequencies in a chosen window. Typically, the numerical handling 
of matrices with approximate dimension 100 x 100 is sufficient even when the 
number of frequencies in the signal C(s) is about 10000 or more. The knowl- 
edge of the operator Ci itself is not required as for a properly chosen basis 
the matrix elements of Cl can be expressed only in terms of the signal C(t). 
The advantage of the filter-diagonalization technique is its numerical stability 
with respect to both the length and complexity (the number and density of 
the contributing frequencies) of the signal. Details of the method of harmonic 
inversion by filter-diagonalization are given in Ref. [70] and in Appendix A.l. 

We now want to apply harmonic inversion as a tool for the high precision anal- 
ysis of quantum spectra. In quantum calculations the bound state spectrum 
is given as a sum of (^-functions, 



(T(w)=J28(w-Wk). (2.21) 



Instead of using the Fourier transform (Eq. 2.7) we want to adjust g qra (w) 
to the functional form of Gutzwiller's semiclassical trace formula (2.3), which 
can be written as 

1 oo 

g sc (w) = g (w) + — J2 E Ko, r e"^ + A^e"^^} . (2.22) 

Z7T po r =l 



The fluctuating part of the semiclassical density of states (2.22) has exactly 
the functional form of the signal C(t) in Eq. 2.14 with t replaced by the 
scaling parameter w. The amplitudes and frequencies {dk,Uk} in Eq. 2.14 
are the amplitudes and scaled actions {A po ^/2ir, — rs po } of the periodic orbit 
contributions and their conjugate pairs {A* r /2TT,rs*}. In order to obtain 
g qm (w) on an evenly spaced grid the spectrum is regularized by convoluting 
it with a narrow Gaussian function having the width a <C l/s max , where s max 
is the scaled action of the longest orbit of interest. The regularized density of 
states reads 

QTH =Ek(f-«fc) = -^Y.^ {W ~ Wkf/2a2 > (2-23) 

k V2na k 
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and is the starting point for the harmonic inversion procedure. The step width 
for the discretization of g^ m (w) is typically chosen as Aw ~ a/3. The convo- 
lution of g qra (w) with a Gaussian function does not effect the frequencies, i.e., 
the values obtained for the scaled actions s po of the periodic orbits, but just 
results in a small damping of the amplitudes 

^po,r ► ^po,r = ^po,rC p ° • (2.24) 



For open systems the density of states is given by 

g ^( w ) = -- Im y , (2.25) 

TV *~£ W -W k 



with complex resonances w k - If the minimum of the resonance widths T^ = 
— 2Im-u;fc is larger than the step width Aw, there is no need to convolute 
£> qm (-u;) with a Gaussian function and the density of states can directly be 
analyzed by harmonic inversion in the same way as for bound systems. 



2. 3 Precision check of the periodic orbit theory 



As a first application of harmonic inversion for the high resolution analysis 
of quantum spectra we investigate the hydrogen atom in a magnetic field 
given by the Hamiltonian (2.11) and compare the results of the harmonic 
inversion method to the conventional Fourier transform presented in Fig. 4. 
We analyze the density of states (2.23) with a = 0.0015 at constant scaled 
energy E = —0.1. From the discussion of the harmonic inversion technique and 
especially Eq. 2.19 it follows that the frequencies in g qm (w), i.e., the actions 
of periodic orbits, are not obtained from a continuous frequency spectrum 
(or a spectrum defined on an equidistant grid of frequencies) but are given 
as discrete and in general complex eigenvalues. Actions with imaginary part 
significantly below zero indicate "ghost" orbit contributions to the quantum 
spectrum. The actions s^/2ti and absolute values of the amplitudes \A qm \ are 
given in the first three columns of Table 1. The last two columns present the 
classical actions Sp /27r of the periodic orbits and the absolute values of the 
semiclassical amplitudes 



„prim In 

| A ci| = s p° IJ1 ) (2.26) 

'det(M po -/)| 

with sP" m the action of the primitive periodic orbit, i.e., the first recurrence of 
the orbit. A graphical comparison between the Fourier transform, the high res- 
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olution quantum recurrence spectrum, and the semiclassical recurrence spec- 
trum is presented in Fig. 5. The crosses in Fig. 5b are the (complex) actions 
obtained by the harmonic inversion of the quantum mechanical density of 
states, and the squares are the actions of the (real) classical periodic orbits. 
The crosses and squares are in excellent agreement with a few exceptions, 
e.g., around s/2ix ps 2.2 and s/2ir ps 2.8 which will be discussed later. The 
amplitudes of the periodic orbit contributions are illustrated in Fig. 5a. The 
solid sticks are the amplitudes obtained by harmonic inversion of the quantum 
spectrum and the dashed sticks (hardly visible under solid sticks) present the 
corresponding semiclassical results. For comparison, the conventional Fourier 
transform recurrence spectrum is drawn as a solid line. To visualize more 
clearly the improvement of the high resolution recurrence spectrum compared 
to the conventional Fourier transform a small part of the recurrence spec- 
trum (Fig. 5a) around s/27r = 2.6 is enlarged in Fig. 6. The smooth line is the 
absolute value of the conventional Fourier transform. Its shape suggests the ex- 
istence of at least three periodic orbits but obviously the recurrence spectrum 
is not completely resolved. The results of the high resolution spectral analy- 
sis are presented as sticks and crosses at the positions defined by the scaled 
actions s^™ with peak heights |A qm |. Note that the positions of the peaks are 
considerably shifted with respect to the maxima of the conventional Fourier 
transform. To compare the quantum recurrence spectrum with Gutzwiller's 
periodic orbit theory the semiclassical results are presented as dashed sticks 
and squares in Fig. 6. For illustration the shapes of periodic orbits are also 
shown (in semiparabolic coordinates ji = \Jr + z, v = \Jr — z) . For these 
three orbits the agreement between the semiclassical and the high resolution 
quantum recurrence spectrum is nearly perfect, deviations are within the stick 
widths. The relative deviations between the quantum mechanical and classical 
actions and amplitudes are given in Table 2. 

The excellent agreement between the classical periodic orbit data and the 
quantum mechanical results obtained by harmonic inversion of the density of 
states may appear to be in contradiction to the fact that Gutzwiller's trace 
formula is an approximation, i.e., only the leading term of the semiclassical 
H expansion of the periodic orbit sum. The small deviations given in Table 2 
are certainly due to numerical limitations and do not indicate effects of higher 
order h corrections. The reason is that the higher order h contributions of the 
periodic orbits do not have the functional form of Eq. 2.22 as a linear super- 
position of exponential functions of the scaling parameter w. As the harmonic 
inversion procedure adjusts the quantum mechanical density of states g qm (w) 
to the ansatz (2.22) with free parameters s po and A po , the exact parameters s po 
and A po of the classical periodic orbits should provide the optimal adjustment 
to the quantum spectrum within the lowest order % approximation. However, 
the harmonic inversion technique can also be used to reveal the higher or- 
der ft, contributions of the periodic orbits in the quantum spectra, as will be 
demonstrated in Section 2.6. 
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The quantum mechanical high resolution recurrence spectrum is not in good 
agreement with the classical calculations around s/2tt pa 2.2 (see crosses 
and squares in Fig. 5b). This period is close to the scaled action s po /2n = 
(—2E)~ l l 2 = 2.2361 of the parallel orbit, which undergoes an infinite series 
of bifurcations with increasing energy [100] . If orbits are too close to bifurca- 
tions they can no longer be treated as isolated periodic orbits, which in turn 
is the basic assumption for Gutzwiller's semiclassical trace formula (2.3). The 
disagreement between the crosses and squares in Fig. 5 around s/2n pa 2.2 
thus indicates the breakdown of the semiclassical ansatz (2.3) for non-isolated 
orbits. The structures around s/2n pa 2.8 can also not be explained by real 
classical periodic orbits. They are related to uniform semiclassical approxima- 
tions and complex ghost orbits as will be discussed in the next section. 



2.4 Ghost orbits and uniform semiclassical approximations 



Gutzwiller's periodic orbit theory (2.3) is valid for isolated periodic orbits 
where the determinant det(Mp — I), i.e., the denominator in Eq. 2.5 is nonzero. 
However, the periodic orbit amplitudes A po ^ r diverge close to the bifurcation 
points of orbits where ML, — / has a vanishing determinant. To remove the 
unphysical singularity from the semiclassical expressions all periodic orbits 
which participate at the bifurcation must be considered simultaneously in 
a uniform semiclassical approximation [40,41]. The uniform solutions can be 
constructed with the help of catastrophe theory [38,39], and have been studied, 
e.g., for the kicked top [36,42,43,101] and the hydrogen atom in a magnetic 
field [37,44]. In the vicinity of bifurcations "ghost" orbits, i.e., periodic orbits 
in the complex continuation of phase space can be very important. In general 
the ghost orbits have real or complex actions, s po . As can be shown from 
the asymptotic expansions of the uniform semiclassical approximations [37,44] 
those ghosts with positive imaginary part of the action, Ims po > 0, are of 
physical relevance. They contribute as 

A iRcs po w _— ImspoW 
^ipo e e 

to Gutzwiller's periodic orbit sum (2.3), i.e., the modulations of the ghost 
orbits are exponentially damped with increasing scaling parameter w. 

Here we will only present a brief derivation of uniform semiclassical approxi- 
mations. For details we refer the reader to the literature [36,37,42-44,101]. Our 
main interest is to demonstrate how bifurcations and ghost orbits can directly 
be uncovered in quantum spectra with the help of the harmonic inversion 
technique and the high resolution recurrence spectra. Note that a detailed 
investigation of these phenomena is in general impossible with the Fourier 
transform because the orbits participating at the bifurcation have nearly the 
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same period and thus cannot be resolved in the conventional recurrence spec- 
tra. We will discuss two different types of bifurcations by way of example of the 
hydrogen atom in a magnetic field: The hyperbolic umbilic and the butterfly 
catastrophe. 

In the following section we will adopt the symbolic code of Ref. [102] for the 
nomenclature of periodic orbits. Introducing scaled semiparabolic coordinates 
jx = 7 1 / 3 ^/r + z and v = ^ 1 / 3 v / r — z the scaled Hamiltonian of the hydrogen 
atom in a magnetic field reads 

h = \(pI+pI) - £(/i 2 + v 2 ) + i(A 2 + A 4 ) = 2 . (2.27) 



As pointed out in [102] the effective potential V{ji, u) = —E{ji 2 + u 2 ) + (/i 4 z/ 2 + 
yU 2 z/ 4 )/8 is bounded for E > by four hyperbolas. When the smooth potential 
is replaced with hard walls periodic orbits can be assigned with the same 
ternary symbolic code as orbits of the four-disk scattering problem (see Fig. 
7). In a first step a ternary alphabet {OAC} is introduced. The symbol labels 
scattering to the opposite disk, and the symbols C and A label scattering to 
the neighboring disk in clockwise or anticlockwise direction, respectively. The 
orbit in Fig. 7 is coded 0C0A. The ternary code can be made more efficient with 
respect to the exchange symmetry of the C and A symbol if it is redefined in 
the following way. The Cs and As are replaced with the sympol +, if consecutive 
letters ignoring the 0s are equal and they are replaced with the symbol -, if 
consecutive letters differ. With this redefinition the orbit shown in Fig. 7 is 
coded 0-0-, or because of its periodicity even simpler 0-. For the hydrogen 
atom in a magnetic field there is a one to one correspondence between the 
periodic orbits and the ternary symbolic code at energies E > E c = +0.3287 
[19,21]. Below the critical energy orbits undergo bifurcations and the code 
is not unique. However, nearly all unstable periodic orbits can be uniquely 
assigned with this code even at negative energies. 



2.4-1 The hyperbolic umbilic catastrophe 

As a first example we investigate the structure at s po /2tt = 2.767 in the recur- 
rence spectrum of the hydrogen atom in a magnetic field at constant scaled 
energy E = E^~ 2 ^ = —0.1 (see Fig. 5 and Table 1). At this energy no pe- 
riodic orbit with an action close to s po /2n = 2.767 does exist. However, the 
strong recurrence peak in Fig. 5a indicates a near bifurcation. This bifurca- 
tion and the corresponding semiclassical approximation have been studied in 
detail in Refs. [44,80]. Four periodic orbits are created through two nearby 
bifurcations around the scaled energy E ~ —0.096 where we search for both 
real and complex "ghost" orbits. For the nomenclature of the real orbits we 
adopt the symbolic code of Ref. [102] as explained above. At scaled energy 
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E^ = —0.09689, the two orbits 00+- and +++ are born in a tangent bi- 
furcation. At energies E < E b , a prebifurcation ghost orbit and its complex 
conjugate exist in the complex continuation of the phase space. Orbit 00+- is 
born unstable, and turns stable at the slightly higher energy E^ = —0.09451. 
This is the bifurcation point of two additional orbits, 0-+ — and its time rever- 
sal +, which also have ghost orbits as predecessors. The graphs of the real 

orbits at energy E = are shown as insets in Fig. 8, and the classical periodic 
orbit parameters are presented as solid lines in Figs. 8 and 9. Fig. 8 shows the 
difference in scaled action between the orbits. The action of orbit 0-+ — (or its 

time reversal +), which is real also for its prebifurcation ghost orbits, has 

been taken as the reference action. The uniform semiclassical approximation 
for the four orbits involved in the bifurcation can be expressed in terms of the 
diffraction integral of a hyperbolic umbilic catastrophe 

+oo +oo 

q?(x,y) = f dp f dqj*^***) ( 2 .28) 



with 

$(p, q; x, y) = p 3 + q 3 + y(p + q) 2 + x{p + q) . (2.29) 



For our convenience the function $(p, q; x, y) slightly differs from the standard 
polynomial of the hyperbolic umbilic catastrophe given in Ref. [39] but the 
diffraction integral (2.28) can be easily transformed to the standard represen- 
tation. The four stationary points of the integral (2.28) are readily obtained 
from the condition V$ = as 



Po = -<?o = ±J-x/3 => $(p , go; x,y) = (2.30) 



and 



2 _, 4 2 X 

4/8 \ /4 x\ 3 / 2 

Hp ,qo;x,y) = -yi^-y 2 -xJT4:i K -y 2 --j . (2.31) 

The function &(po, qo; x, y) must now be adapted to the classical action of the 
four periodic orbits, i.e., AS = wAs ps $(po, qo] x, y), which is well fulfilled for 



x = aw 2 ' 3 



(E-Ei 2) ) ; y = bw 1 l\ (2.32) 



27 



and constants a = —5.415, b = 0.09665, as can be seen from the dashed lines 
in Fig. 8. Note that the agreement holds for both the real and complex ghost 
orbits. 

The next step to obtain the uniform solution is to calculate the diffraction 
integral (2.28) within the stationary phase approximation. For E > E b there 
are four real stationary points (po,<7o) (see Eqs. 2.30 and 2.31), and after 
expanding <&(p,q;x,y) around the stationary points up to second order in p 
and q, the diffraction integral becomes the sum of Fresnel integrals, viz. 



,., r ;n 2n ^ 7re 



!v§» 2 -*=F4 f^-ir ±f 



9(x, y) ^ -p= + E "T =— . (2.33) 

V-3a; + _ 7(4^2 _ 3a .) T 2y^Ay 2 - 3x 



The terms of Eq. 2.33 can now be compared to the standard periodic orbit 
contributions (2.26) of Gutzwiller's trace formula. In our example the first 
term is related to the orbit 0-+ — (with a multiplicity factor of 2 for its time 

reversal +), and the other two terms are related to the orbits 00+- and 

+++ for the upper and lower sign, respectively. The phase shift in the 

numerators describe the differences of the action AS and of the Maslov index 
A/x = =pl relative to the reference orbit 0-+ — . The denominators are, up to 
a factor cw 1 ^ 3 , with c = 0.1034, the square root of | det(M — I)\, with M the 
stability matrix. Fig. 9 presents the comparison for the determinants obtained 
from classical periodic orbit calculations (solid lines) and from Eqs. 2.32 and 
2.33 (dashed lines). The agreement is very good for both the real and complex 
ghost orbits, similar to the agreement found for As in Fig. 8. The constant 
c introduced above determines the normalization of the uniform semiclassical 
approximation for the hyperbolic umbilic bifurcation which is finally obtained 
as [44,80] 

A unilorm (E,w) = {c/-K)s w l l^(aw 2 l\E - E^bw 1 ' 3 ) j^-W , 

(2.34) 

with So and /xo denoting the orbital action and Maslov index of the reference 
orbit 0-+ — , and the constants a, b, and c as given above. 

The comparison between the amplitudes (2.26) of the conventional semiclas- 
sical trace formula for isolated periodic orbits and the uniform approximation 
(2.34) for the hyperbolic umbilic catastrophe is presented in Fig. 10 at the 
magnetic field strengths 7 = 10~ 7 , 7 = 10~ 8 , and 7 = 10~ 9 . For graphical 
purposes we suppress the highly oscillatory part resulting from the function 
exp[i(S /H — f/U )] by plotting the absolute value of A(E,w) instead of the 
real part. The dashed line in Fig. 10 is the superposition of the isolated peri- 
odic orbit contributions from the four orbits involved in the bifurcations. The 
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modulations of the amplitude are caused by the constructive and destructive 
interference of the real orbits at energies E > E^ and are most pronounced at 
low magnetic field strength (see Fig. 10c). The amplitude diverges at the two 
bifurcation points. For the calculation of the uniform approximation (2.34) 
we numerically evaluated the catastrophe diffraction integral (2.28) using a 
more simple and direct technique as described in [103]. Details of our method 
which is based on Taylor series expansions are given in Appendix CI. The 
solid line in Fig. 10 is the uniform approximation (2.34). It does not diverge 
at the bifurcation points but decreases exponentially at energies E < E b . 
At these energies no real orbits exist and the amplitude in the standard for- 
mulation would be zero when only real orbits are considered. However, the 
exponential tail of the uniform approximation (2.34) is well reproduced by a 
ghost orbit [37,36] with positive imaginary part of the complex action. As can 
be shown, the asymptotic expansion of the diffraction integral (2.28) has, for 
x ^> 0, exactly the form of Eq. 2.26 but with complex action s and determi- 
nant det(M — I). The ghost orbit contribution is shown as dash-dotted line 
in Fig. 10. 

To verify the hyperbolic umbilic catastrophe in the quantum spectrum we 
analyze all three, the exact quantum spectrum, the uniform semiclassical ap- 
proximation (2.34), and Gutzwiller's periodic orbit formula for isolated pe- 
riodic orbits by means of the harmonic inversion technique at scaled energy 
E = —0.1, which is slightly below the two bifurcation points. The part of the 
complex action plane which is of interest for the hyperbolic umbilic catastro- 
phe is presented in Fig. 11. The two solid peaks mark the positions Sk/2ir and 
the absolute values of amplitudes \Ay\ obtained from the quantum spectrum. 
As mentioned above, at this energy only one classical ghost orbit is of physi- 
cal relevance and marked as dash-dotted peak in Fig. 11. The position of that 
peak is in good agreement with the quantum result but the amplitude is en- 
hanced. This enhancement is expected for isolated periodic orbit contributions 
near bifurcations which become singular exactly at the bifurcation points. The 
harmonic inversion analysis of the uniform approximation (2.34) at constant 
scaled energy E = —0.1 in the same range < w < 140 is presented as dashed 
peaks in Fig. 11. The two peaks agree well with the quantum results for both 
the complex actions and amplitudes. The enhancement of the ghost orbit peak 
and the additional non-classical peak observed in the quantum spectrum are 
therefore clearly identified as artifacts of the bifurcation, i.e., the hyperbolic 
umbilic catastrophe. 



2.4-2 The butterfly catastrophe: Uncovering the "hidden" ghost orbits 

We now investigate the butterfly catastrophe, which is of importance, e.g., in 
photoabsorption spectra of the hydrogen atom in a magnetic field. In contrast 
to the density of states the photoabsorption spectra for dipole transitions from 
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a low-lying initial state to highly excited final states can be measured experi- 
mentally [45-47]. The semiclassical photoabsorption cross-section is obtained 
by closed orbit theory [90-93] as the superposition of a smooth background 
and sinusoidal modulations induced by closed orbits starting at and returning 
back to the nucleus. Although the derivation of the semiclassical oscillator 
strength for dipole transitions (closed orbit theory) differs from the derivation 
of the semiclassical density of states (periodic orbit theory) the final results 
have a very similar structure and therefore spectra can be analyzed in the 
same way by conventional Fourier transform or the high resolution harmonic 
inversion technique. In closed orbit theory the semiclassical oscillator strength 
is given by 

f = f° + f° sc (2.35) 

with /° the oscillator strength of the field free hydrogen atom at energy E = 
and 



9(271-^3/2 

r c = 2(E f -E i )w- l ' 2 \m ]T ^yU^ 

ci.o. k y \rh12 I 

x y m ($i,k)y m ($f,k) expli (s k w + -mf k - -fi° k + -J I (2.36) 



the fluctuating part of the oscillator strength. In (2.36) Ef and E^ are the 
energies of the final and initial state, m is the magnetic quantum number, 
Sfe, Tj, = 7Tfc, and /i° are the scaled action, scaled time, and the Maslov 
index of closed orbit k, i?j and i?/ are the starting and returning angle of the 
orbit with respect to the magnetic field axis, and rhyi is an element of the 
scaled monodromy matrix. The angular functions y m ($) depend on the initial 
state and polarization of light (see Appendix B). For more details see Refs. 
[47,91,93]. The fluctuating terms (2.36) of the photoabsorption cross section 
are sinusoidal functions of the scaling parameter w = 7" 1 / 3 despite a factor of 
w -i/2 Tq obtain the same functional form as Eq. 2.22 for the density of states, 
which is required for harmonic inversion, we multiply the oscillator strength / 
by w x l 2 for both the semiclassical and quantum mechanical photoabsorption 
spectra. 

In analogy to Gutzwiller's trace formula Eq. 2.36 for photoabsorption spectra 
is valid only for isolated orbits and diverges at bifurcations of orbits, where 
rh\2 is zero. Near bifurcations the closed orbit contributions in (2.36) must be 
replaced with uniform semiclassical approximations, which have been studied 
in detail for the fold, cusp, and butterfly catastrophe in Refs. [37,81]. Here we 
restrict the discussion to the butterfly catastrophe, which is especially interest- 
ing because of the existence of a "hidden" ghost orbit, which can be uncovered 
in the photoabsorption spectrum of the hydrogen atom in a magnetic field with 



30 



the help of the harmonic inversion technique. As an example we investigate 
real and ghost orbits related to the period doubling of the perpendicular or- 
bit R\. (For the closed orbits we adopt the nomenclature of Refs. [45,46], 
see also Fig. 2.) This closed orbit bifurcation is more complicated because 
various orbits with similar periods undergo two different elementary types of 
bifurcations at nearly the same energy. The structure of bifurcations and the 
appearance of ghost orbits can clearly be seen in the energy dependence of 
the starting angles i?j in Fig. 12a. Two orbits R\ and R\^ are born in a saddle 

-0.31735345, ^ = 1.3465. Below the bifurcation 
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(1) 



node bifurcation at E, 

energy we find an associated ghost orbit and its complex conjugate. Orbit R 
is real only in a very short energy interval (AE ps 0.001), and is then involved 
in the next bifurcation at E^ — —0.31618537, i?j = tt/2. This is the period 
doubling bifurcation of the perpendicular orbit R±, which exists at all energies 
(i?i = 7r/2 in Fig. 12a). The real orbit Rl b separates from R\ at energies below 
the bifurcation point, i.e. a real orbit vanishes with increasing energy. Conse- 
quently, associated ghost orbits are expected at energies above the bifurcation, 
i.e. E > El , and indeed such "postbifurcation" ghosts have been found. Its 
complex starting angles are also shown in Fig. 12a. The energy dependence of 
scaled actions, or, more precisely, the difference As with respect to the action 
of the period doubled perpendicular orbit R 2 , is presented in Fig. 12b (solid 
lines), and the graph for the monodromy matrix element m 12 is given in Fig. 
12c. It can be seen that the actions and the monodromy matrix elements of 
the ghost orbits related to the saddle node bifurcation of R\ become complex 
at E < E^ , while these parameters remain real for the postbifurcation ghosts 
at E > El . The two bifurcations are so closely adjacent that formulae for the 
saddle node bifurcation and the period doubling are no reasonable approxima- 
tion to As(E) and fhi2(E) in the neighborhood of the bifurcations. However, 
both functions can be fitted well by the more complicated formulae [37] 



As = k 



a 



E-E, 



(2) 



!«± 



a(E-M 2) ) + l 



3/2 



(2.37) 



and 



m 12 



M E-E, 



(2) 



m 12 = 4:M(E-E. 



i(2) 



(orbit i?2 

AM 



a 



l±Jo-[E-Ei 2) )-\ 



(2.38) 



(Rl,R}, b , and ghosts) 



with k = 3.768 • 10~ 4 , a = 763.6, and M = 13.52 (see dashed lines in Figs. 
12b and 12c). Note that Eqs. (2.37) and (2.38) describe the complete scenario 
for the real and the ghost orbits including both the saddle node and period 
doubling bifurcations. We also mention that orbits with angles di ^ df have 
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to be counted twice because they correspond to different orbits when traversed 
in either direction, and therefore a total number of five closed orbits, includ- 
ing ghosts, is considered here in the bifurcation scenario around the period 
doubling of the perpendicular orbit. 

The bunch of trajectories forming the butterfly is given by the Hamilton- Jacobi 
equations (with p^ = dS/d^ and p v = dS/drf) 



3dp 5 v + icp\p% + p\ \a - to)p v - m] = ° 

where £ and rj are rotated semiparabolic coordinates 



(2.39) 



£ = n cos -f + i/sm-T- 

■&t . tit 

r\ — v cos -p- — fi sin -+ 



2r cos 



■&--&J 



2r sin — - — - 



(2.40) 
(2.41) 



so that the £ and i] axes are now parallel and perpendicular to the returning 
orbit. The parameters c, d, and £o i n Eq. 2.39 will be specified later. With 
VnJVi ~ aAf9j we obtain 



r){i) = 3a 5 d(A$ t ) 5 + 2a 3 c(A^) 3 + a(A^)(£ - £ ) • 



(2.42) 



The butterfly is illustrated in Fig. 13. Depending on the number of real solu- 
tions Ai?j of Eq. (2.42) there exist one, three, or five orbits returning to each 
point (£,r)). The different regions are separated by caustics. 

The uniform semiclassical approximations for closed orbits near bifurcations 
can in general be written as 



f° k sc = 2(E - E t ) y/antii sm^ f y m (^)y m ^ f ) 
xlmj^lexp (i '" v " ■ ■''' 
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where the complex amplitude A is defined implicitly by the diffraction integral 

[37,81] 
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(2.44) 



To find a uniform semiclassical approximation for the butterfly catastrophe 
we have to solve the Hamilton- Jacobi equations (2.39) at least in the vicinity 
of the central returning orbit and to insert the action S(r, i?) into Eq. (2.44). 
For the classical action we obtain 



AS\r, 0) = ±V& + ±$,(0 - ^) 2 + i c (0 - 7 ) 4 - lrf(0 - d f f , 

(2.45) 

and for the determinant in the denominator of (2.44) we find ±a\/8r in the 
limit r ^> and i? ~ 0/. Summing up in (2.44) the contributions of the 
incoming and the outgoing orbit (with Maslov indices A/i A as for the cusp) 
we obtain the integral I(r) and the amplitude A for the butterfly catastrophe 



I(r) = 2 1 /4 r -i/4| a |~i/2 exp (_ -^ cog J^ _ 3 \ 
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where 
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fit 



(2.48) 



is an analytic function in both variables x and y. Its numerical calculation and 
asymptotic properties are discussed in Appendix C.2. The uniform result for 
the oscillatory part of the transition strength now reads 



r c = 2{E _ Ei) v / sin ^ sin ^ / 3; m (^)3; m ^ / )2 5 /2 7r |a|- 1 / 2 rf- 1 / 6 



xlm< exp i 






n c 

2" 



* (-d- 1/3 qo , -cd~ 



2/3 



(2.49) 



It is very illustrative to study the asymptotic behavior of the uniform ap- 
proximation (2.49) as we obtain, on the one hand, the relation between the 
parameters a, c, d, and £o and the actions and the monodromy matrix ele- 
ments of closed classical orbits, and, on the other hand, the role of complex 
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ghost orbits related to this type of bifurcation is revealed. In the following 
we discuss both limits £ 3> 0, i.e. scaled energy E 3> E^ , and £ <C 0, i.e. 



E^E, 



(i) 



(2) 



Asymptotic behavior at scaled energy E ~^> E\ ] Applying Eq. (CIO) 
from Appendix C.2 to the \l/-function in the uniform approximation (2.49), 
we obtain the asymptotic formula for £o ^ 
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(2.50) 



Comparing with the solutions for isolated closed orbits (2.36) we can identify 
the contributions of three real closed orbits. The classical action of orbit 1 is 
S^j, its Maslov index is /i° and the monodromy matrix element mu is given 
by 



(i) 



-a& = - 7 - 1 / 3 M( J B - E b 



(2.51) 



where the parameter M can be determined by closed orbit calculations (see 
Eqs. 2.37 and 2.38). Orbits 2 and 3 are symmetric with respect to the z = 
plane and have the same orbital parameters, i.e. Maslov index /x° + 1 and the 
classical action and the monodromy matrix element 
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(2.53) 



~ (2) 

In the example of the bifurcation of orbit R 2 at scaled energy El = —0.31618537, 
orbit 1 is the orbit R 2 perpendicular to the magnetic field axis, while orbits 2 
and 3 can be identified with R\ traversed in both directions (i? i;3 = it — i? i)2 )- 
With the help of Eqs. (2.51) to (2.53) and classical scaling properties of the ac- 
tion and the monodromy matrix the parameters a, c, d, and £o in the uniform 
approximation (2.49) can now be expressed completely in terms of closed-orbit 
parameters k, a, and M (see Eqs. 2.37 and 2.38), 
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and the uniform approximation for the butterfly catastrophe finally reads 
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In the classical analysis complex ghost orbits were discovered both below and 
above the bifurcation energy. At E > E^ they have the property that the 
classical action and the monodromy matrix remain real, although coordinates 
and momenta in phase space are complex. These ghost orbits do not appear 
in the asymptotic expansion (2.50) of the uniform approximation (2.49), and 
therefore, in analogy with the cusp catastrophe (see [37]), they do not have 
a physical meaning. The situation is different at energy E < E b where a 
"hidden ghost" with physical meaning will be revealed in the following. 



Asymptotic behavior and "hidden ghost" at scaled energy E <C E[ 

At scaled energies below the bifurcation point we can apply the asymptotic 
formula (C.12) from Appendix C.2 to the \l/-function in the uniform approxi- 
mation (2.49) and obtain 
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The first term in Eq. (2.56) can be identified as a real orbit with the same 
classical action and monodromy matrix element as in (2.50), but with a Maslov 
index increased by one. The second term in (2.56) is a ghost orbit contribution 
resulting from a superposition of two closed orbits with complex action and 
monodromy matrix element, 
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traversed in both directions. The positive imaginary part of the classical ac- 
tion results in an exponential damping of the ghost orbit contribution to the 
oscillator strength amplitude with decreasing energy similar to the situation 
at a fold catastrophe. In contrast to the fold catastrophe (see Ref. [37]) and to 
the hyperbolic umbilic catastrophe discussed in Section 2.4.1 the ghost orbit 
related to a butterfly catastrophe is always accompanied by a real orbit with 
almost the same classical action. Because the contribution of the real orbit 
is not exponentially damped its amplitude at energies where the asymptotic 
formula (2.56) is valid is much stronger than the ghost contribution. Therefore 
we call the ghost orbit in (2.56) a "hidden ghost" [37,81]. Note that classi- 
cally the complex conjugate of the hidden ghost orbit also exists. The negative 
imaginary part of its classical action would result in an unphysical exponential 
increase of amplitude with decreasing energy, and consequently the complex 
conjugate ghost orbit does not appear in the asymptotic formula (2.56). In 
Ref. [37] it was admitted that it might be rather difficult to find evidence for 
hidden ghost orbits in the Fourier transform of experimental or theoretical 
scaled energy spectra. However, with the harmonic inversion technique we are 
now able to uncover the hidden ghosts in high resolution quantum recurrence 
spectra. 



Uncovering the hidden ghost orbit To uncover the hidden ghost orbit 
related to the period doubling bifurcation of the perpendicular orbit in pho- 
toabsorption spectra of the hydrogen atom in a magnetic field we calculated 
the quantum spectrum at constant scaled energy E = —0.35, which is suffi- 
ciently far below the two bifurcation energies around E = —0.317 so that the 
real orbit R2 and the prebifurcation ghost orbit of R\ are approximately iso- 
lated orbits. We calculated 2823 transitions from the initial state |2p0) to final 
states with magnetic quantum number m = in the region w = 7 -1 / 3 < 100. 
The scaled photoabsorption spectrum was analyzed by conventional Fourier 
transform and by the high resolution harmonic inversion technique. The in- 
teresting part of the recurrence spectrum in the region 1.86 < s/2n < 1.93 
is presented in Fig. 14. The conventional Fourier transform (smooth line) has 
a maximum at s/2tt = 1.895 which is roughly twice the period of the per- 
pendicular orbit but does not give any hint on the existence of a ghost orbit. 
The key points of the harmonic inversion analysis are that the resolution of 
the recurrence spectrum is not restricted by the uncertainty principle of the 
Fourier transform and that the method supplies complex frequencies s qm of 
the analyzed quantum spectrum / qm (-u;), which can be interpreted as complex 
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actions Sp of ghost orbits. The high resolution spectral analysis uncovers one 
real and two complex actions s qm around s/2ir = 1.9, which are marked as 
crosses in Fig. 14b. They can be compared to actions of real and complex 
classical orbits. We find two closed orbits, i.e., the period doubling of the per- 
pendicular orbit with scaled action s^ /2tt = 1.896011 and a ghost orbit with 
scaled action Sp /27r = 1.894401 — iO. 006372, which are marked as squares in 
Fig. 14b. The shapes of the real and complex closed orbits are presented as in- 
sets in Fig. 14a (in semiparabolic coordinates //, v). The real part of the ghost 
orbit (solid line) is similar to the shape of orbit R\ which is created as a real 
orbit at much higher energy E = —0.317. The actions and also the amplitudes 
in the quantum and classical recurrence spectrum (see crosses and squares in 
Fig. 14) agree very well for the two closed orbit recurrences. The deviation 
between the imaginary parts of the actions of the quantum and classical ghost 
orbit can be explained by the fact that the orbits are only approximately iso- 
lated at energy E = —0.35. The correct semiclassical formula is the uniform 
semiclassical approximation (2.55) of the butterfly catastrophe. This interpre- 
tation is supported by the occurrence of a third complex resonance with large 
imaginary part of the action and small amplitude in the quantum recurrence 
spectrum (Fig. 14), which has no classical analogue. The non-classical peak 
is similar to the occurrence of a non-classical peak in Fig. 11 for the hyper- 
bolic umbilic catastrophe and is clearly an artifact of the near bifurcation. 
The uncovering of the hidden ghost orbit, which obviously is impossible with 
the conventional Fourier transform (solid line in Fig. 14) demonstrates that 
harmonic inversion is a very powerful tool for the analysis of quantum spectra, 
which can reveal structures and information from the spectra that have been 
unattainable before. 



2.5 Symmetry breaking 



Methods for high resolution spectral analysis are also helpful for a direct ob- 
servation of symmetry breaking effects in quantum spectra. In many cases 
physical systems possess symmetries, e.g., a cylindrical symmetry around a 
fixed axis. In such situations the closed or periodic orbits of the classical sys- 
tem appear in continuous families. All members of the family have the same 
stability parameters and periods, i.e., they are observed as one peak in the 
Fourier transform recurrence spectrum of the quantum mechanical density of 
states or the transition spectrum. The dynamics is profoundly changed when 
the symmetry is broken by a (weak) external perturbation. The behavior of 
the dynamics and the corresponding semiclassical theory is described in Ref. 
[104]. In general, out of a continuous family of orbits only two closed or peri- 
odic orbits survive the symmetry breaking. Recently, quantum manifestations 
of symmetry breaking have been observed experimentally for atoms in external 
fields: The cylindrical symmetry of the hydrogen atom in a magnetic field was 
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broken in crossed magnetic and electric fields [54], and a "temporal symmetry 
breaking" was studied on lithium atoms in an oscillating electric field [55,56]. 

However, the expected splitting of recurrence peaks cannot be observed di- 
rectly because of the finite resolution of the Fourier transform. This can be 
seen in Fig. 15 for the hydrogen atom in crossed fields. Fig. 15a presents 
segments of the experimental recurrence spectra at constant scaled energy 
E = —0.15 and scaled electric field strengths < / = F7" 4 / 3 < 0.055. F 
is the electric field strength in atomic units (Fq = 5.14 x 10 9 V/cm). For 
comparison Fig. 15b shows the theoretical recurrence spectra obtained from 
semiclassical closed orbit theory [54]. Spectra have been Fourier transformed 
in the range 34.0 < w = 7 -1 / 3 < 61.7. The recurrence structure corresponds to 
the classical orbit in the plane perpendicular to the magnetic field axis drawn 
in Fig. 16. At vanishing electric field trajectories starting at the origin are ex- 
actly closed at each return to the nucleus (thin line in Fig. 16a). For nonzero 
electric field, only two orbits with slightly different periods (shown with heavy 
lines in Fig. 16a) return exactly to the origin. A closeup of the returning part 
of the trajectories as they approach the nucleus at scaled energy E = —0.15 
and field strength / = 0.012 is given in Fig. 16b. The two closed orbits return 
to the nucleus diagonally; all others follow near-parabolic paths near the nu- 
cleus. In Fig. 15 no splitting of the recurrence peak can be observed due to 
the finite resolution of the Fourier transform spectra. However, the amplitude 
of the peak changes with increasing electric field strength, and indeed in Ref. 
[54] the symmetry breaking was identified indirectly by the constructive and 
destructive interference of the two orbits resulting in a Bessel function type 
modulation of amplitudes of recurrence peaks as a function of the strength of 
the symmetry breaking perturbation. 

We now want to apply harmonic inversion to directly uncover the splitting 
of recurrence peaks when symmetries are broken. As in Ref. [54] we investi- 
gate the hydrogen atom in crossed magnetic and electric fields. Schrodinger's 
equation for hydrogen in crossed fields reads 



2P 2 - - + lf L z + g7 V + y 2 ) + Fx 



tf = Em , (2.59) 



with 7 and F the magnetic and electric field strength (in atomic units), and 
E the energy. The eigenvalue problem (2.59) was solved numerically for fixed 
external fields 7 and F in Refs. [96,105]. To study the effects of symmetry 
breaking on the quantum spectra we want to use the scaling properties of the 
classical system and to analyze spectra at constant scaled energy E = E^~ 2 ^ 
and scaled electric field / = F^~ A ^. Eigenvalues are obtained for the scaling 
parameter w = 7~ 1 / 3 . Introducing dilated coordinates f = 7 2 / 3 r Eq. 2.59 reads 



w 2 



P 2 



* + w [L z ] * + 



hx 2 + f) + 2fx-2E-- 
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* = . (2.60) 



The peculiarity of (2.60) is that it is a quadratic (instead of linear) eigenvalue 
equation for the scaling parameter w, which cannot be solved straightforwardly 
with standard diagonalization routines for linear eigenvalue problems. To solve 
the quadratic Schrodinger equation (2.60) we use the following technique. We 
write 



w 



A = A + AA 



(2.61) 



with An fixed and assume that AA is small. Now the scaling parameter w in 
the paramagnetic term of (2.60) can be approximated by 



w 



Vx 



A + 



1 



2V\ 



-.AX. 



(2.62) 



Replacing w in (2.60) with this approximation yields a generalized linear eigen- 
value problem for AA, 



AoP 2 + \fx~oL z + \{x 2 + f) + 2fx-2E-- 
v 4 r 
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L ? 



*, 



(2.63) 



which can be solved numerically with the Spectral Transformation Lanczos 
Method (STLM) [106]. The numerical details for the diagonalization of the 
linearized Schrodinger equation (2.63) in a complete set of Sturmian type 
basis functions are similar to the diagonalization of Eq. 2.59 at constant field 
strengths as described in [96]. It is convenient to choose identical parameters 
for both An and the center of the spectral transformation, so that the STLM 
method provides eigenvalues in the local neighborhood of A ~ An, i.e., AA is 
small. However, the obtained eigenvalues are still not very precise because of 
the approximation (2.62). In a second step we therefore apply perturbation 
theory and use the eigenvectors |\P) of Eq. 2.63 to solve the quadratic equation 
(2.60) separately for each eigenvalue, w n , in the corresponding one-dimensional 
subspace, i.e. 



w 2 n (#„|p 2 |^„) + w n (y n \L z \V n ) 



1 

4 ( 



+ (^ n \^-(x 2 + y 2 ) 



2fx-2E--\^ n ) = 0. 



(2.64) 



By this procedure the accuracy of the eigenvalues is significantly improved. An 
alternative method for the exact solution of the quadratic eigenvalue equation 
(2.60) is described in Ref. [107]. 
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We now investigate the symmetry breaking in the hydrogen atom at constant 
scaled energy E = —0.5. Without electric field (/ = 0) the two shortest 
closed orbit recurrences are the perpendicular orbit at s/2tt = 0.872 and the 
parallel orbit at s/2-k = 1. When the cylindrical symmetry is broken by a 
crossed electric field, all three-dimensional orbits should split into two nearby 
peaks. The only exception is the orbit parallel to the magnetic field axis, 
which does not appear as a continuous family of closed orbits. To verify the 
symmetry breaking in quantum spectra we calculated the photoabsorption 
spectrum (transitions from the initial state \2p0) to the final states with even 
z-parity) up to w = 50 for the crossed field atom at constant scaled energy E = 
—0.5 and scaled field strength / = 0.02. The interesting part of the resulting 
recurrence spectrum obtained by both the conventional Fourier transform and 
the high resolution harmonic inversion technique are presented in Fig. 17. 
The conventional Fourier transform (smooth line) shows two peaks around 
s/27r = 0.87 (the perpendicular orbit) and s/2-k = 1 (the parallel orbit). 
However, none of the peaks is split. The high resolution recurrence spectrum 
obtained by harmonic inversion is drawn as solid sticks and crosses in Fig. 17 
and clearly exhibits a splitting of the recurrence peak of the perpendicular 
orbit at s/2ix m 0.87. For comparison the semiclassical recurrence spectrum 
is presented as dashed sticks and squares. In the plane perpendicular to the 
magnetic field axis two closed orbits have been found with slightly different 
classical actions s cl /2n = 0.864 and s cl /2n = 0.881. The shapes of these orbits 
are illustrated as insets in Fig. 17. As can be seen the semiclassical and the high 
resolution quantum recurrence spectrum are in excellent agreement, i.e., the 
symmetry breaking in crossed fields has been directly uncovered by harmonic 
inversion of the quantum mechanical photoabsorption spectrum. 



2.6 H expansion of the periodic orbit sum 



In the previous sections we have introduced harmonic inversion as a powerful 
method for the high resolution analysis of quantum spectra, which allows a 
direct and quantitative comparison of the quantum spectra with semiclassical 
theories. However, the excellent agreement to many significant digits for both 
the periods and amplitudes of quantum mechanical recurrence peaks with the 
semiclassical periodic orbit contributions (see Section 2.3) may be surprising 
for the following reason. Periodic orbit theory is exact only for a special class 
of systems, e.g., the geodesic motion on a surface with constant negative cur- 
vature [8]. In general, the semiclassical periodic orbit sum is only the leading 
order contribution of an infinite series in powers of the Planck constant [22- 
24]. Therefore it might be expected that only the high resolution analysis of 
semiclassical spectra yields perfect agreement with periodic orbit theory while 
the harmonic inversion of quantum spectra should show small but noticeable 
deviations from the semiclassical recurrence spectra. As will be shown in the 
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following, the absence of such deviations is related to the functional form of the 
h expansion of the periodic orbit sum and special properties of the harmonic 
inversion method. 

In scaling systems, where the classical action of periodic orbits scales as 

Sp /h = Spo/^cflf = S po W , 

the scaling parameter plays the role of an effective Planck constant, i.e. 

w = Ki , 

and the h expansion of the periodic orbit sum can therefore be written as a 
power series in w" 1 . The fluctuating part of the semiclassical response function 
is given by 



g[w) 



OO OO 1 

E 9n(w) = E — E A^e is ^ . (2.65) 



n=0 n=0 W po 



The A!j$ are the complex amplitudes of the n th order periodic orbit contri- 
butions including phase information from the Maslov indices. When quantum 
spectra are analyzed in the semiclassical regime, i.e., at sufficiently high scal- 
ing parameter w, the higher order h correction terms (n > 1) of the series 
(2.65) are certainly small compared to the zeroth order terms. However, the 
reason why the higher order contributions are not uncovered by the harmonic 
inversion method is that only the zeroth order (n = 0) terms in (2.65) fulfill 
the ansatz (2.14) required for harmonic inversion, i.e., they are exponential 
functions with constant frequencies and amplitudes. The higher order terms 
(n > 1) have amplitudes decreasing ~ w~ n with increasing scaling parame- 
ter, and thus do not fulfill the ansatz (2.14). Therefore only the zeroth order 
amplitudes A p °] can be obtained as converged parameters from the high res- 
olution harmonic inversion analysis of the quantum spectra. The terms with 
n > 1 have similar properties as weak "noise" [70] and are separated by the 
harmonic inversion method from the "true" signal. In other words, the zeroth 
order approximation of the periodic orbit sum (2.65) with amplitudes A p ^ 
given by Gutzwiller's trace formula is the best fit to the quantum spectra 
within the given ansatz as a linear superposition of exponential functions of 
the scaling parameter w, and the corresponding parameters are obtained from 
the harmonic inversion procedure. 

However, the higher order terms of the % expansion (2.65) can be revealed 
by harmonic inversion as will be demonstrated in the following. The periodic 
orbit terms A p ^ can be obtained provided that the quantum spectrum and 
the (n — l) st order eigenvalues Wk n -i are given. We can then calculate the 
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difference between the exact quantum mechanical and the (n — l) st order 
response function 



n— 1 oo oo i 



S qm M - E 9jM = E 9iM = E ~ 3 T,^y Sp ° W ■ (2-66) 

j=0 j=n j=n po 



The leading order terms in Eq. 2.66 are ~ w n , i.e., multiplication with w n 
yields 



w n 



n-1 

^ qm W-EftW 

i=o 



J2A$e iSpoW + o(-} . (2.67) 

po \ W / 



In Eq. 2.67 we have restored the functional form (2.14), i.e., a linear superpo- 
sition of exponential functions of w. The harmonic inversion of the function 
(2.67) will now provide the periods s po and the n th order amplitudes A!j$ of the 
h expansion (2.65). We will illustrate the method on two different examples, 
i.e., the circle billiard and the three disk scattering problem. 



2.6.1 The circle billiard 

The circle billiard is an integrable and even separable bound system, and 
has been chosen here mainly for the sake of simplicity, since all the relevant 
physical quantities, i.e., the quantum and semiclassical eigenenergies, and the 
periodic orbits can easily be obtained. In polar coordinates (p, 0) and after 
separation of the 0-motion Schrodinger's equation reads 

h 2 (\ d d m 2 \ 

-pi- ~ -r ^m{p) = Em m {p) , (2.68) 



2M Xpdp^dp p 2 



with M the mass of the particle and m the angular momentum quantum 
number. The wave functions must fulfill the boundary condition ty m (R) = 
with R the radius of the circle billiard. Defining 

h 2 k 2 

E = 

2M 



and after substitution of z = kp Schrodinger's equation (2.68) is transformed 
into the differential equation for the Bessel functions [108] 

z 2 J^(z) + zJ' m (z) + (z 2 - m 2 )J m (z) = . (2.69) 
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The quantum mechanical eigenvalues k c ^ n are obtained from the boundary 
condition J m (/c^it!) = as zeros of the Bessel functions. We calculated 
numerically the first 32469 eigenenergies of the circle billiard in the region 
k n ^ m R < 510. Note that states with m/0 are twofold degenerate. 

The semiclassical eigenenergies are obtained from an EBK torus quantiza- 
tion or, after separation of the 0-motion, even more simply from the one- 
dimensional WKB quantization of the radial motion in the centrifugal poten- 
tial, i.e. 

R R 

3 s 
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2 / p p dp = 


--2hk f 


m/k 


m/k 



y/l - {m/kpfdp = 2irh (n + jj , (2.70) 



with n = 0,1,2,... the radial quantum number. The r.h.s. of Eq. 2.70 takes 
into account the correct boundary conditions of the semiclassical wave func- 
tions at the classical turning points. The zeroth order semiclassical eigenener- 
gies are finally obtained from the quantization condition 

I 777- f e \\ 

kR\j\ — [m/kR) 2 — \m\ arccos — — - = n ( n -\ — 1 . (2-71) 



The semiclassical spectrum has been calculated with the help of Eq. 2.71 in 
the same region k^ m R < 510 as the exact quantum spectrum. 

For the comparison of the spectra with periodic orbit theory we need to calcu- 
late the periodic orbits and their physical quantities. For the circle billiard all 
quantities are obtained analytically. In the following we choose R — l. The pe- 
riodic orbits of the circle billiard are those orbits for which the angle between 
two bounces is a rational multiple of 2ir, i.e., the periods £ po are obtained from 
the condition 



po 



2m r sin7, (2.72) 



with 



7 = 7T 



rritj, = 1,2,... the number of turns of the orbit around the origin, and m r = 
2171^,2171^ + 1, . . . the number of reflections at the boundary of the circle. 
Periodic orbits with m r ^ 2vn^ can be traversed in two directions and thus 
have multiplicity 2. Because the classical dynamics of the circle billiard is 
regular the Berry- Tabor formula [25] must be applied instead of Gutzwiller's 
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trace formula for the calculation of the semiclassical density of states. We 
obtain 

-I 1 OO OO 1 

-7r9{k) = -/= E 9n(k) = E Vn E^^ , (2-73) 

VK V fv n=0 n= o ^ po 



which basically differs from Eq. 2.65 by a factor of A;" 1 / 2 on the l.h.s. of (2.73). 
[Note that for billiard systems the scaling parameter is the absolute value of 
the wave vector, w = k = \p\/h, and the action is proportional to the length of 
the orbit, S po = hk£ po .} The zeroth order periodic orbit amplitudes obtained 
from the Berry- Tabor formula read 

4o } = fM^^ ) > ( 2 - 74 ) 



with /i po = 3m r the Maslov index. For the calculation of the first order periodic 
orbit amplitudes A$ in Eq. 2.73 we adopt the method of Alonso and Gaspard 
[23] . After a lengthy calculation we finally obtain 

^ = IV^r 5 " 2 3 S ;f 7 e "" fMP °" f } • ( 2 - 75 ) 

I 3sm ' 7 

/(7) 



A detailed derivation of Eqs. 2.74 and 2.75 will be given elsewhere [87,88]. 

When the quantum mechanical density of states, or, more precisely, the spec- 
trum k~ 1 ' 2 g qm (k) is analyzed by the harmonic inversion method the periodic 
orbit quantities £ po and ,4^ are obtained to very high precision [87]. However, 
we are now interested in the % expansion of the periodic orbit sum and want 
to verify the first order corrections A„ 1 } directly in the quantum spectrum. We 
therefore analyze the difference spectrum Ag(k) = g qra (k) — g sc (k) between 
the quantum and the semiclassical density of states. A small part of this spec- 
trum at k ~ 200 is presented in Fig. 18. The absolute values of the peak 
heights mark the multiplicities of the states. To restore the functional form 
which is required for the harmonic inversion procedure the difference spec- 
trum was multiplied by \fk, and the resulting signal VkAg(k) was analyzed 
in the region 100 < k < 500. The results are presented in Fig. 19. For a direct 
comparison with Eq. 2.75 we transform for each periodic orbit the obtained 
period £ po and the first order amplitude ,4^ into the quantities 

7 = arcsin 2lt and /(7)s vi^ | 4 S| - 
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These quantities are plotted as crosses in Fig. 19. The periodic orbits are 
marked by the numbers (m^,,m r ). The solid line 

m 5 ~ 2sm2 ~> 



3 sin 3/2 



7 



is the result of Eq. 2.75. As can be seen the theoretical curve and the crosses 
obtained by harmonic inversion of the difference between the quantum and 
semiclassical density of states are in excellent agreement. 

The semiclassical accuracy of the zeroth order eigenenergies of the circle bil- 
liard has been discussed in Refs. [109,110] in terms of, e.g., the average error in 
units of the mean level spacing. The analysis presented here provides a more 
physical interpretation of the deviations between the quantum and semiclas- 
sical eigenenergies in terms of the higher order corrections of the h expanded 
periodic orbit sum (2.65). The direct application of the series (2.65) for the 
semiclassical quantization beyond the Gutzwiller and Berry- Tabor approxi- 
mation will be discussed in Section 3.6. 



2.6.2 The three disk scattering system 

As a second example for the investigation of the H expansion of the periodic 
orbit sum we now consider a billiard system consisting of three identical hard 
disks with unit radius R—l, displaced from each other by the same distance, 
d. The classical, semiclassical, and quantum dynamics of this scattering system 
has been studied by Gaspard and Rice [111]. In recent years the system has 
served as a prototype model for periodic orbit quantization by cycle expansion 
techniques [9-11,14,15]. If the disks are separated by a distance d > 2.0481419 
there is a one-to-one identity between the periodic orbits and a symbolic code, 
whereas for d < 2.04821419 pruning of orbits sets in [112]. The geometry of 
the three disk scattering system is shown in Fig. 20. The symbolic code of 
a periodic orbit is obtained by numbering of the disks, '1', '2', and '3', and 
by bookkeeping the reflections of the orbit at the three disks. E.g., 121313232 
is the symbolic code of the primitive orbit with cycle length n = 9 drawn 
with a dashed line in Fig. 20. The bar, which is often omitted, indicates the 
periodicity of the orbit. The three disk scattering system is invariant under 
the symmetry operations of the group Cs v , i.e., three reflections at symmetry 
lines and two rotations by 27r/3 and 47r/3. Periodic orbits have a v symmetry if 
they are invariant under reflections and C3, C\ symmetry if they are invariant 
under rotations. After symmetry decomposition the periodic orbits in the fun- 
damental domain (see Fig. 20b) can be classified by a binary symbolic code of 
symbols '0' and '1', where each '0' represents a change between clockwise and 
anticlockwise scattering in the original three disk system [9]. The symbolic 
code of the orbit 121313232 in Fig. 20a restricted to the fundamental domain 



45 



is 100100100 or because of its periodicity even simpler 100. The quantum 
resonances are also classified by symmetries. Resonances with Ai (A 2 ) sym- 
metry are symmetric (antisymmetric) under reflections at a symmetry line, 
and resonances with E symmetry are invariant (up to a complex phase factor 
c with c 3 = 1) under rotations by 27r/3 and An/3. The resonances in the E 
subspace are twofold degenerate. In the following we analyze resonances with 
A\ symmetry. 

The three disk scattering system does not have any bound states. However, 
the fluctuating part of the density of states can be written as 

e(*) = ~ Im Err (2 - 76) 



with k n the complex resonances (poles) of the Green's function. For distance 
d — 6 the quantum mechanical A\ resonances have been calculated by Wirzba 
[14,15,113] in the region < Re k < 250, and are plotted with the + symbols in 
Fig. 21c. The corresponding fluctuating part of the quantum density of states 
is shown in Fig. 21a. This spectrum can be analyzed by harmonic inversion 
to extract the lengths (periods) of the orbits and the zeroth order amplitudes. 
However, as in the previous section for the circle billiard we here want to 
investigate the difference between the quantum and semiclassical spectrum. 
The semiclassical Ai resonances have been calculated by Wirzba using the 
12 th order in the curvature expansion of the Gutzwiller-Voros zeta-function 
[14,15,113] and are plotted by the crosses (x) in Fig. 21c. The differences 
between the quantum and semiclassical resonances are usually much smaller 
than the size of the symbols in Fig. 21c and are visible only at small values of 
Re k. However, the deviations become clearly pronounced in Fig. 21b, which 
presents the difference between the quantum and zeroth order semiclassical 
density of states 

k[g qm (k) - g sc (k)] = -ImV Me il »° k + O ( r) . (2.77) 



Eq. 2.77 implies that by harmonic inversion of the spectrum in Fig. 21b we can 
extract the first order amplitudes AQ} of the % expansion of the periodic orbit 
sum. In Table 3 we present the ratio of the zeroth and first order amplitudes, 
l«4poV*4po | hl obtained by harmonic inversion of the spectra in Fig. 21a and 21b 
for all periodic orbits up to cycle length 5 in the symmetry reduced symbolic 
code [9]. These quantities have been calculated by Alonso and Gaspard [23] 
and we present their results |*4poA4po 1°' in Table 3 for comparison. In Ref. 
[23] the periodic orbits have not been symmetry reduced and the set of orbits 
chosen by Alonso and Gaspard is complete only up to cycle length 3 in the 
symmetry reduced symbolic code. However, Table 3 clearly illustrates a very 
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good agreement between the amplitudes obtained by the harmonic inversion 
analysis of spectra and the results of Ref. 



We have here obtained information about the first order terms of the % expan- 
sion of the periodic orbit sum from the analysis of the quantum and the lowest 
(zeroth) order semiclassical resonances. In general, the quantum spectrum and 
its n th order semiclassical approximation are required to extract information 
about the in + l) st order terms of the % expansion series (2.65). 



3 Periodic orbit quantization by harmonic inversion 



Since the development of periodic orbit theory by Gutzwiller [7,8] it has be- 
come a fundamental question as to how individual semiclassical eigenenergies 
and resonances can be obtained from periodic orbit quantization for classi- 
cally chaotic systems. A major problem is the exponential proliferation of the 
number of periodic orbits with increasing period, resulting in a divergence of 
Gutzwiller's trace formula at real energies and below the real axis, where the 
poles of the Green's function are located. The periodic orbit sum is a Dirichlet 
series 



JW = E^*> t 3 - 1 ) 



where the parameters A n and s n are the amplitudes and periods (actions) 
of the periodic orbit contributions. In most applications Eq. 3.1 is absolutely 
convergent only in the region Im w > c > with c the entropy barrier of 
the system, while the poles of g(w), i.e., the bound states and resonances, are 
located on and below the real axis, Im w < 0. Thus, to extract individual 
eigenstates, the semiclassical trace formula (3.1) has to be analytically con- 
tinued to the region of the quantum poles. Up to now no general procedure is 
known for the analytic continuation of a non-convergent Dirichlet series of the 
type of Eq. 3.1. All existing techniques are restricted to special situations. For 
bound and ergodic systems the semiclassical eigenenergies can be extracted 
with the help of a functional equation and the mean staircase function (Weyl 
term), resulting in a Riemann-Siegel look-alike formula [59-61]. Alternative 
semiclassical quantization conditions based on a semiclassical representation 
of the spectral staircase [65,66] and derived from a quantum version of a clas- 
sical Poincare map [63] are also restricted to bound and ergodic systems. For 
systems with a symbolic dynamics the periodic orbit sum (3.1) can be re- 
formulated as an infinite Euler product, which can be expanded in terms of 
the cycle length of the symbolic code. If the contributions of longer orbits are 
shadowed by the contributions of short orbits the cycle expansion technique 
can remarkably improve the convergence properties of the series and allows 
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to extract the bound states and resonances of bound and open systems, re- 
spectively [9-11,58]. A combination of the cycle expansion technique with a 
functional equation for bound systems has been studied by Tanner et al. [114]. 
However, the existence of a simple symbolic code is restricted to very few sys- 
tems, and cycle expansion techniques cannot be applied, e.g., to the general 
class of systems with mixed regular-chaotic classical dynamics. 

In this chapter we present a general technique for the analytic continuation 
and the extraction of poles of a non-convergent series of the type of Eq. 3.1. 
The method is based on harmonic inversion by filter-diagonalization. The ad- 
vantage of the method is that it does not depend on special properties of the 
system such as ergodicity or the existence of a symbolic dynamics for periodic 
orbits. It does not even require the knowledge of the mean staircase function, 
i.e., the Weyl term in dynamical systems. The only assumption we have to 
make is that the analytic continuation of the Dirichlet series g(w) (Eq. 3.1) 
is a linear combination of poles (w — u>fc) _1 , which is exactly the functional 
form of, e.g., a quantum mechanical response function with real and complex 
parameters Wk representing the bound states and resonances of the system, 
respectively. To demonstrate the general applicability and accuracy of our 
method we will apply it to three systems with completely different proper- 
ties, first the zeros of the Riemann zeta function [115,116], as a mathematical 
model for a bound system, second the three disk scattering system as a phys- 
ical example of an open system with classically chaotic dynamics, and third 
the circle billiard as an integrable system. 

As pointed out by Berry [117] the density of zeros of Riemann's zeta function 
can be written, in formal analogy with Gutzwiller's semiclassical trace formula, 
as a non- convergent series, where the "periodic orbits" are the prime numbers. 
A special property of this system is the existence of a functional equation 
which allows the calculation of Riemann zeros via the Riemann- Siegel for- 
mula [115-118]. An analogous functional equation for quantum systems with 
an underlying chaotic (ergodic) classical dynamics has served as the basis for 
the development of a semiclassical quantization rule for bound ergodic systems 
[59-61]. The Riemann zeta function has also served as a mathematical model 
to study the statistical properties of level distributions [118-120]. We will 
demonstrate in Section 3.1 that harmonic inversion can reveal the Riemann 
zeros with extremely high accuracy and with just prime numbers as input 
data. The most important advantage of our method is, however, its wide ap- 
plicability, i.e., it can be generalized in a straightforward way to non-ergodic 
bound or open systems, and the procedure for periodic orbit quantization by 
harmonic inversion will be discussed in Section 3.2. 

As an example of periodic orbit quantization of a physical system we inves- 
tigate in Section 3.3 the three disk scattering problem, which is an open and 
non-ergodic system. Its classical dynamics is purely hyperbolic, and the peri- 
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odic orbits can be classified by a complete binary symbolic code. This system 
has served as the prototype for the development of cycle expansion techniques 
[9-1 1] . When applying the harmonic inversion technique to the three disk scat- 
tering system we will highlight the general applicability of our method by not 
having to make use of its symbolic dynamics in any way. 

The power of the method will be illustrated in Section 3.4 on a challenging 
physical system which has not been solved previously with any semiclassi- 
cal quantization technique. The hydrogen atom in a magnetic field shows a 
transition from near-integrable to chaotic classical dynamics with increasing 
excitation energy. We apply periodic orbit quantization by harmonic inversion 
to the hydrogen atom in the mixed regular-chaotic regime. 

An important question is the efficiency of methods for periodic orbit quan- 
tization, i.e., the number of periodic orbits required for the calculation of a 
certain number of poles of the response function g(w). It is evident that meth- 
ods invoking special properties of a given system may be remarkably efficient. 
E.g., the Riemann-Siegel type formulae [59-61] require the periodic orbits up 
to a maximum period which is by about a factor of four shorter compared to 
the required signal length for the harmonic inversion technique discussed in 
Section 3.2. We therefore propose in Section 3.5 an extension of the harmonic 
inversion method to the harmonic inversion of cross- correlated periodic orbit 
sums. The method uses additional semiclassical information obtained from 
a set of linearly independent smooth observables and allows to significantly 
reduce the number of periodic orbits and thus to improve the efficiency of 
periodic orbit quantization by harmonic inversion. 

Periodic orbit theory yields exact eigenenergies only in exceptional cases, e.g., 
for the geodesic motion on the constant negative curvature surface [8]. As 
already discussed in Section 2.6, Gutzwiller's periodic orbit sum is, in gen- 
eral, just the leading order term of an infinite series in powers of the Planck 
constant. Methods for the calculation of the higher order periodic orbit con- 
tributions were developed in [22-24]. In Section 3.6 we demonstrate how the 
higher order % corrections of the periodic orbit sum can be used to improve 
the accuracy of the semiclassical eigenenergies, i.e., to obtain eigenenergies 
beyond the standard semiclassical approximation of periodic orbit theory. 

Both methods, the harmonic inversion of cross-correlated periodic orbit sums, 
which allows to significantly reduce the required number of periodic orbits 
for semiclassical quantization, and the calculation of eigenenergies beyond 
the lowest order H approximation will be illustrated in Section 3.7 by way of 
example of the circle billiard. Furthermore, for this system we will calculate 
semiclassically the diagonal matrix elements of various operators. The circle 
billiard is an integrable system and therefore the Berry- Tabor formula for 
integrable systems [25,26] is valid in this case rather than Gutzwiller's trace 
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formula. However, periodic orbit quantization by harmonic inversion can be 
applied for both regular and chaotic systems as well which demonstrates the 
universality and wide applicability of the method. 

Finally, in Section 3.8 we will calculate semiclassically the photoabsorption 
spectra of atoms in external fields. Applying a combination of closed orbit 
theory [90-93], the cross-correlation approach of Section 3.5, and the harmonic 
inversion method we obtain individual non-diagonal transition matrix elements 
between low-lying initial states and strongly perturbed Rydberg states of the 
magnetized hydrogen atom. 



3.1 A mathematical model: Riemann's zeta function 



Our goal is to introduce our method for periodic orbit quantization by har- 
monic inversion using, as an example, the well defined problem of calculating 
zeros of the Riemann zeta function. There are essentially two advantages of 
studying the zeta function instead of a "real" physical bound system. First, 
the Riemann analogue of Gutzwiller's trace formula is exact, as is the case 
for systems with constant negative curvature [8,66], whereas the semiclassical 
trace formula for systems with plane geometry is correct only to first order 
in h. This allows a direct check on the precision of the method. Second, no 
extensive periodic orbit search is necessary for the calculation of Riemann 
zeros, as the only input data are just prime numbers. It is not our intention 
to introduce yet another method for computing Riemann zeros, which, as an 
objective in its own right, can be accomplished more efficiently by specific pro- 
cedures. Rather, in our context the Riemann zeta function serves primarily as 
a mathematical model to illustrate the power of our technique when applied 
to bound systems. 



3.1.1 General remarks 

Before discussing the harmonic inversion method we start with recapitulating 
a few brief remarks on Riemann's zeta function necessary for our purposes. 
The hypothesis of Riemann is that all the non-trivial zeros of the analytic 
continuation of the function 

oo _, 

(( z ) = J2 n~ z = J] (l - P~ z ) ~ , (Re z > 1, p : primes) (3.2) 

71=1 P 

have real part |, so that the values w = Wk, defined by 

C Q - iw^j = 0, (3.3) 
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are all real or purely imaginary [115,116]. The Riemann staircase function for 
the zeros along the line z — | — iw, defined as 

oo 

N(w) = '£G(w-w k ), (3.4) 

fe=i 



i.e. the number N(w) of zeros with w^ < w, can be split [115-117] into a 
smooth part, 



N(w) = — arg T ( — | — iw ) In n + 1 

7r \4 2 / 2tt 



and a fluctuating part, 



■^osc(w) = lim Imln( -- i(w + irj) ) . (3.6) 

7r ^;^o \2 



Substituting the product formula (3.2) (assuming that it can be used when 
Re z — |) into (3.6) and expanding the logarithms yields 

1 °° 1 

^oscM = -- Im £ £ — ^ e<™W . (3.7) 



Therefore the density of zeros along the line z — \ — iw can formally be written 
as 

Qosc{w) = —j^- = lmg(w) (3.8) 

aw 7T 



with the response function g(w) given by the series 

sM = «ES^ mln(p) > ( 3 - 9 ) 



V m=\ 



P" 



which converges only for Imw > |. Obviously Eq. 3.9 is of the same type 
as the response function (3.1), with the entropy barrier c = |, i.e., Eq. 3.9 
does not converge on the real axis, where the Riemann zeros are located. The 
mathematical analogy between the above equation and Gutzwiller's periodic 
orbit sum 

fo sc (£)^--rm]T„4 po e^ , (3.10) 

*T po 
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with A po the amplitudes and S po the classical actions (including phase in- 
formation) of the periodic orbit contributions, was already pointed out by 
Berry [117,59]. For the Riemann zeta function the primitive periodic orbits 
have to be identified with the primes p, and the integer m formally counts the 
"repetitions" of orbits. The "amplitudes" and "actions" are then given by 



ln(p) 

p Ti 
S pm = mw\n(p) . (3-12) 



A -i-^- (3 11 

■V - ' m/2 ' ^° 



Both equation (3.9) for the Riemann zeros and - for most classically chaotic 
physical systems - the periodic orbit sum (3.10) do not converge. In partic- 
ular, zeros of the zeta function, or semiclassical eigenstates, cannot be ob- 
tained directly using these expressions. The problem is to find the analytic 
continuation of these equations to the region where the Riemann zeros or, for 
physical systems, the eigenenergies and resonances, are located. Eq. 3.9 is the 
starting point for our introduction and discussion of the harmonic inversion 
technique for the example of the Riemann zeta function. The generalization 
of the method to periodic orbit quantization (Eq. 3.10) in Section 3.2 will be 
straightforward. 

Although Eq. 3.9 is the starting point for the harmonic inversion method, 
for completeness we quote the Riemann- Siegel formula, which is the most 
efficient approach to computing Riemann zeros. For the Riemann zeta function 
it follows from a functional equation [115] that the function 



.1 1 \ 1 ■ 

argr — I — iw wmn 

1 4 2 J 2 



C[\~iw) (3.13) 



Z{w) = exp <^ —i 
is real, and even for real w. The asymptotic representation of Z(w) for large 



lnt[y/w/2n] f JCt/ \ 11 

ry, \ v-> cos{7riv (it; J — w In n) 

Z[w) = -2 ^ -— 

71=1 n 

_ ( _ 1) i^V^]( 2 _lY cos( 27 r(t 2 -t-l/16)) 

v ; \wj cos(27rf) '" ' K ' J 

with t = Jw/2tx — Int [J w/2tc] is known as the Riemann-Siegel formula and 
has been employed (with several more correction terms) in effective methods 
for computing Riemann zeros [118]. Note that the principal sum in (3.14) 



has discontinuities at integer positions of Jw/2ir, and therefore the Riemann 
zeros obtained from the principal sum are correct only to about 1 to 15 per- 
cent of the mean spacing between the zeros. The higher order corrections to 
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the principal Riemann-Siegel sum remove, one by one, the discontinuities in 
successive derivatives of Z(w) at the truncation points and are thus essential 
to obtaining accurate numerical results. An alternative method for improving 
the asymptotic representation of Z{w) by smoothing the cutoffs with an error 
function and adding higher order correction terms is presented in [61]. An 
analogue of the functional equation for bound and ergodic dynamical systems 
has been used as the starting point to develop a "rule for quantizing chaos" 
via a "Riemann-Siegel look-alike formula" [59-61]. This method is very effi- 
cient as it requires the least number of periodic orbits, but unfortunately it 
is restricted to ergodic systems on principle reasons, and cannot be general- 
ized either to systems with regular or mixed classical dynamics or to open 
systems. By contrast, the method of harmonic inversion does not have these 
restrictions. We will demonstrate that Riemann zeros can be obtained directly 
from the "ingredients" of the non-convergent response function (3.9), i.e., the 
set of values A pm and S pm , thus avoiding the use of the functional equation, the 
Riemann-Siegel formula, the mean staircase function (3.5), or any other spe- 
cial property of the zeta function. The comparison of results in Section 3.1.3 
will show that the accuracy of our method goes far beyond the Riemann-Siegel 
formula (3.14) without higher order correction terms. The main goal of this 
section is to demonstrate that because of the formal equivalence between Eqs. 
(3.9) and (3.10) our method can then be applied to periodic orbit quantization 
of dynamical systems in Section 3.2 without any modification. 



3.1.2 The ansatz for the Riemann zeros 

To find the analytic continuation of Eq. (3.9) in the region Im w < | we 
essentially wish to fit g(w) to its exact functional form, 

V w - w k + iO 



arising from the definition of the Riemann staircase (3.4). The "multiplicities" 
dk in Eq. 3.15 are formally fitting parameters, which here should all be equal 
to 1. It is hard to directly adjust the non- convergent (on the real axis) series 
g(w) to the form of g cx (w). The first step towards the solution of the problem 
is to carry out the adjustment for the Fourier components of the response 
function, 



1 r ■ °° Info) 

C(s) = —J g(w)e—dw = i£ £ -^ S(s - mln(p)) , (3.16) 



which after certain regularizations (see below) is a well-behaved function of 
s. Due to the formal analogy with the results of periodic orbit theory (see 
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Eqs. 3.11 and 3.12), C(s) can be interpreted as the recurrence function for the 
Riemann zeta function, with the recurrence positions s pm = S pm /w = m ln(p) 
and recurrence strengths of periodic orbit returns A pm = i\n(p)p~ m ^ 2 . The 
exact functional form which now should be used to adjust C(s) is given by 

i + r 

Cex(s) = — J g cx (w)e~ lsw dw = -i £ d k e~^ s . (3.17) 



_ x k— 1 



C cx (s) is a superposition of sinusoidal functions with frequencies Wk given by 
the Riemann zeros and amplitudes d k — 1. (It is convenient to use the word 
"frequencies" for w k referring to the sinusoidal form of C(s). We will also use 
the word "poles" in the context of the response function g(w).) 

Fitting a signal C(s) to the functional form of Eq. 3.17 with, in general, both 
complex frequencies Wk and amplitudes dk is known as harmonic inversion, 
and has already been introduced in Section 2.2 for the high resolution analysis 
of quantum spectra. The harmonic inversion analysis is especially non-trivial 
if the number of frequencies in the signal C(s) is large, e.g., more than a 
thousand. It is additionally complicated by the fact that the conventional way 
to perform the spectral analysis by studying the Fourier spectrum of C(s) will 
bring us back to analyzing the non- convergent response function g(w) defined 
in Eq. 3.9. Until recently the known techniques of spectral analysis [67] would 
not be applicable in the present case, and it is the filter-diagonalization method 
[68-70] which has turned the harmonic inversion concept into a general and 
powerful computational tool. 

The signal C(s) as defined by Eq. 3.16 is not yet suitable for the spectral 
analysis. The next step is to regularize C(s) by convoluting it with a Gaussian 
function to obtain the smoothed signal, 



+oo 

C a {s) = -L^- f C(s')e-( s - s ' )2/2(T2 d S ' 
V2^a J v ' 

— oo 



= * VV 1 ^ e -(- m ln( P ))W (3 . 18) 

that has to be adjusted to the functional form of the corresponding convolution 
of C ex (s). The latter is readily obtained by substituting dk in Eq. 3.17 by the 
damped amplitudes, 

d k ^d k a) =d k e- w " a2/2 . (3.19) 



The regularization (3.18) can also be interpreted as a cut of an infinite num- 
ber of high frequencies in the signal which is of fundamental importance for 
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numerically stable harmonic inversion. Note that the convolution with the 
Gaussian function is no approximation, and the obtained frequencies Wk and 
amplitudes dk corrected by Eq. 3.19 are still exact, i.e., do not depend on o. 
The convolution is therefore not related to the Gaussian smoothing devised 
for Riemann zeros in [121] and for quantum mechanics in [122], which provides 
low resolution spectra only. 

The next step is to analyze the signal (3.18) by harmonic inversion. The con- 
cept of harmonic inversion by filter-diagonalization has already been explained 
in Section 2.2 and the technical details are given in Appendix A.l. Note that 
even though the derivation of Eq. 3.18 assumed that the zeros Wk are on the 
real axis, the analytic properties of C a (s) imply that its representation by Eq. 
3.18 includes not only the non-trivial real zeros, but also all the trivial ones, 
Wk = —i(2k + |), k — 1, 2, . . ., which are purely imaginary. The general har- 
monic inversion procedure does not require the frequencies to be real. Both 
the real and imaginary zeros Wk will be obtained as the eigenvalues of the 
non-Hermitian generalized eigenvalue problem, Eq. A. 7 in Appendix A.l. 



3.1.3 Numerical results 

For a numerical demonstration we construct the signal C a (s) using Eq. 3.18 in 
the region s < ln(10 6 ) = 13.82 from the first 78498 prime numbers and with a 
Gaussian smoothing width o = 0.0003. Parts of the signal are presented in Fig. 
22. Up to s ~ 8 the Gaussian approximations to the (^-functions do essentially 
not overlap (see Fig. 22a) whereas for s ^> 8 the mean spacing As between 
successive (^-functions becomes much less than the Gaussian width a = 0.0003 
and the signal fluctuates around the mean C(s) = ie s ^ 2 (see Fig. 22b). From 
this signal we were able to calculate about 2600 Riemann zeros to at least 
12 digit precision. For the small generalized eigenvalue problem (see Eq. A. 7 
in Appendix A.l) we used matrices with dimension J < 100. Some Riemann 
zeros Wk, the corresponding amplitudes dk, and the estimated errors e (see 
Eq. A. 12 in Appendix A.l) are given in Tables 4 and 5. The pole of the zeta 
function yields the smooth background C(s) = ie s ^ 2 of the signal C(s) (see 
the dashed line in Fig. 22). Within the numerical error the Riemann zeros are 
real and the amplitudes are consistent with dk — 1 for non-degenerate zeros. 
To fully appreciate the accuracy of our harmonic inversion technique we note 
that zeros obtained from the principal sum of the Riemann- Siegel formula 
(3.14) deviate by about 1 to 15 percent of the mean spacing from the exact 
zeros. Including the first correction term in (3.14) the approximations to the 
first five zeros read w x = 14.137, w 2 = 21.024, w 3 = 25.018, w 4 = 30.428, and 
w§ = 32.933, which still significantly deviates from the exact values (see Table 
4). Considering even higher order correction terms the results will certainly 
converge to the exact zeros. However, the generalization of such higher order 
corrections to ergodic dynamical systems is a nontrivial task and requires, 
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e.g., the knowledge of the terms in the Weyl series, i.e., the mean staircase 
function after the constant [61,123]. The perfect agreement of our results for 
the Wk with the exact Riemann zeros to full numerical precision is remarkable 
and clearly demonstrates that harmonic inversion by filter-diagonalization is 
a very powerful and accurate technique for the analytic continuation and the 
extraction of poles of a non-convergent series such as Eq. 3.1. 

A few Wk have been obtained (see Table 6) which are definitely not located on 
the real axis. Except for the first at w — ijl they can be identified with the 
trivial real zeros of the zeta function at z — — 2n; n — 1,2,... In contrast to 
the nontrivial zeros with real w k , the numerical accuracy for the trivial zeros 
decreases rapidly with increasing n. The trivial zeros w n = —i(2n + |) are 
the analogue of resonances in open physical systems with widths increasing 
with n. The fact that the trivial Riemann zeros are obtained emphasizes the 
universality of our method and demonstrates that periodic orbit quantization 
by harmonic inversion can be applied not only to closed but to open systems 
as well. The decrease of the numerical accuracy for very broad resonances is 
a natural numerical consequence of the harmonic inversion procedure [69,70]. 
The value w = ijl in Table 6 is special because in this case the amplitude is 
negative, i.e., dk = — 1. Writing the zeta function in the form [116] 

S(\-iw) = C Y[(w - w k ) d "A(w, w k ) (3.20) 

1 k 



where C is a constant and A a regularizing function which ensures convergence 
of the product, integer values dk are the multiplicities of zeros. Therefore it 
is reasonable to relate negative integer values with the multiplicities of poles. 
In fact, ((z) has a simple pole at z — \ — iw — 1 consistent with w = ijl in 
Table 6. 



3.1.4 Required signal length 

We have calculated Riemann zeros by harmonic inversion of the signal C a (s) 
(Eq. 3.18) which uses prime numbers as input. The question arises what are 
the requirements on the signal C ff (s), in particular what is the required signal 
length. In other words, how many Riemann zeros (or semiclassical eigenener- 
gies) can be converged for a given set of prime numbers (or periodic orbits, 
respectively)? The answer can be directly obtained from the requirements on 
the harmonic inversion technique. In general, the required signal length s max 
for harmonic inversion is related to the average density of frequencies g(w) by 
[70] 

Smax ~ 4:iig(w) , (3.21) 
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i.e., s max is about two times the Heisenberg length 

SHcisenbcrg = 2tIq{w) . (3.22) 



From Eq. 3.21 the required number of primes (or periodic orbits) can be 
directly estimated as {# primes p | hip < s max } or {# periodic orbits | s po < 
Smax}- For the special example of the Riemann zeta function the required 
number of primes to have a given number of Riemann zeros converged can 
be estimated analytically. With the average density of Riemann zeros derived 
from (3.5), 

g( w )= =—\ n ( — ) 3.23 

^ v ; dw 2n \2n ' v ; 



we obtain 



w \ ( w 

ln(p ma x) = 2 In ( — J => p max = (7^) • (3.24) 



The number of primes with p < p max can be estimated from the prime number 
theorem 

( \ ^ax (w/2tt) 2 

7T Pmax) ~ T~, T = w , v • (3.25) 

lnQOmax) 2 ln(w/27r) 



On the other hand the number of Riemann zeros as a function of w is given 
by Eq. (3.5). The estimated number of Riemann zeros which can be obtained 
by harmonic inversion from a given set of primes is presented in Fig. 23. 
For example, about 80 zeros (w < 200) can be extracted from the short 
signal C ff (s) with s max = ln(1000) = 6.91 (168 prime numbers) in agreement 
with the estimates given above. Obviously, in the special case of the Riemann 
zeta function the efficiency of our method cannot compete with that of the 
Riemann- Siegel formula method (3.14) where the number of terms is given 



by n max = Int [Jw/2n] and, e.g., 5 terms in Eq. 3.14 would be sufficient to 
calculate good approximations to the Riemann zeros in the region w < 200. 
Our primary intention is to introduce harmonic inversion by way of example 
of the zeros of the Riemann zeta function as a universal tool for periodic 
orbit quantization, and not to use it as an alternative method for solving the 
problem of finding most efficiently zeros of the Riemann zeta function. For the 
semiclassical quantization of bound and ergodic systems a functional equation 
can be invoked to derive a Riemann- Siegel look-alike quantization condition 
[59-61]. In this case the required number of periodic orbits can be estimated 
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from the condition 



Smax ~ HQ[W) — Sjjcisenberg i yo.Zb) 



which differs by a factor of 4 from the required signal length (3.21) for har- 
monic inversion. Obviously, the limitation (3.21) is unfavorable for periodic 
orbit quantization because of the exponential proliferation of the number of 
orbits in chaotic systems, and methods for a significant reduction of the signal 
length s max are highly desirable. In fact, this can be achieved by an extension 
of the harmonic inversion technique to cross-correlation signals. We will return 
to this problem in Section 3.5. 



3.2 Periodic orbit quantization 



As mentioned in Section 3.1 the basic equation (3.9) used for the calculation 
of Riemann zeros has the same mathematical form as Gutzwiller's semiclas- 
sical trace formula. Both series, Eq. 3.9 and the periodic orbit sum (3.10), 
suffer from similar convergence problems in that they are absolutely conver- 
gent only in the complex half-plane outside the region where the Riemann 
zeros, or quantum eigenvalues, respectively, are located. As a consequence, 
in a direct summation of periodic orbit contributions smoothing techniques 
must be applied resulting in low resolution spectra for the density of states 
[78,122]. To extract individual eigenstates the semiclassical trace formula has 
to be analytically continued to the region of the quantum poles, and this was 
the subject of intense research during recent years. 

For strongly chaotic bound systems a semiclassical quantization condition 
based on a semiclassical representation of the spectral staircase 

Af(E)=Y,®(E-E n ) (3.27) 



has been developed by Aurich et al. in Refs. [65,66]. They suggest to replace the 
spectral staircase Af(E) with a smooth semiclassical approximation Af sc (E), 
and to evaluate the semiclassical eigenenergies from the quantization condition 



cos {ttKc(E)} =0 . (3.28) 



For the even parity states of the hyperbola billiard the exact spectral staircase 
Af + (E) and the semiclassical approximation JV^(E) evaluated from 101265 
periodic orbits are shown in Fig. 24. The function cos {irAf+(E)} obtained 
from the smooth spectral staircase J\f+(E) in Fig. 24 is presented in Fig. 25. 



According to the quantization condition (3.28) the zeros of this function are 
the semiclassical energies, and are in good agreement with the true quantum 
mechanical energies marked by triangles in Fig. 25. The quantization condition 
(3.28) was also successfully applied to the motion of a particle on various 
Riemann surfaces with constant negative curvature, e.g., Artin's billiard or 
the Hadamard-Gutzwiller model [65,66]. Another technique for bound and 
ergodic systems is to apply an approximate functional equation and generalize 
the Riemann-Siegel formula (3.14) to dynamical zeta functions [59-61]. The 
Riemann- Siegel look-alike formula has been applied, e.g., for the semiclassical 
quantization of the hyperbola billiard [123]. 

These quantization techniques cannot be applied to open systems. However, 
if a symbolic dynamics for the system exists, i.e., if the periodic orbits can 
be classified with the help of a complete symbolic code, the dynamical zeta 
function, given as an infinite Euler product over entries from classical periodic 
orbits, can be expanded in terms of the cycle length of the orbits [9,58]. The 
cycle expansion series is rapidly convergent if the contributions of long orbits 
are approximately shadowed by contributions of short orbits. The cycle ex- 
pansion technique has been applied, e.g., to the three disk scattering system 
[9-11] (see also Section 3.3), the three body Coulomb system [6,124], and to 
the hydrogen atom in a magnetic field [20]. 

It turns out that the cycle expansion of dynamical zeta functions converges 
very slowly for bound systems. Therefore, a combination of the cycle expan- 
sion method with a functional equation has been developed by Tanner et al. 
[114,125]. Applying a functional equation to the dynamical zeta function they 
conclude that the zeros of the real expression 

D(E) = e~ lnM(E) Z(E) + e iw ^ E) Z*(E) , (3.29) 

with M{E) the mean spectral staircase and Z(E) the cycle expanded dynam- 
ical zeta function should be semiclassical approximations to the eigenvalues. 
The semiclassical functional determinant D(E) is illustrated in Fig. 26 for the 
anisotropic Kepler problem [8] and in Fig. 27 for the closed three disk billiard, 
i.e., the system shown in Fig. 20 with touching disks. The vertical bars in Figs. 
26 and 27 mark the exact eigenvalues. 

The existence of a complete symbolic dynamics is more the exception than the 
rule, and therefore cycle expansion techniques cannot be applied, in particu- 
lar, for systems with mixed regular-chaotic classical dynamics. In this section 
we apply the same technique that we used for the calculation of Riemann 
zeros to the calculation of semiclassical eigenenergies and resonances of phys- 
ical systems by harmonic inversion of Gutzwiller's periodic orbit sum for the 
propagator. The method only requires the knowledge of all orbits up to a 
sufficiently long but finite period and does neither rely on an approximate 
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semiclassical functional equation, nor on the existence of a symbolic code for 
the orbits. The universality of the method will therefore allow the investiga- 
tion of a large variety of systems with an underlying chaotic, regular, or even 
mixed classical dynamics. The derivation of an expression for the recurrence 
function to be harmonically inverted is analogous to that in Section 3.1.2. 



3.2.1 Semiclassical density of states 

Following Gutzwiller [7,8] the semiclassical response function for chaotic sys- 
tems is given by 

g sc (E)=gl c (E) + 52A po e iS »°, (3.30) 

po 



where gQ C (E) is a smooth function and the ,Sp and A po are the classical ac- 
tions and weights (including phase information given by the Maslov index) of 
periodic orbit contributions. Eq. (3.30) is also valid for integrable [25,26] and 
near-integrable [27,28] systems but with different expressions for the ampli- 
tudes A po - It should also be possible to include complex "ghost" orbits [36,37] 
and uniform semiclassical approximations [40,44] close to bifurcations of pe- 
riodic orbits in the semiclassical response function (3.30). The eigenenergies 
and resonances are the poles of the response function but, unfortunately, its 
semiclassical approximation (3.30) does not converge in the region of the poles, 
whence the problem is the analytic continuation of g sc (E) to this region. 

In the following we make the (weak) assumption that the classical system has 
a scaling property (see Section 2.1), i.e., the shape of periodic orbits does not 
depend on the scaling parameter, w, and the classical action scales as 

S P o = ws po . (3.31) 



Examples of scaling systems are billiards [9,16], Hamiltonians with homoge- 
neous potentials [126,127], Coulomb systems [6], or the hydrogen atom in 
external magnetic and electric fields [47,20]. Eq. 3.31 can even be applied for 
non-scaling, e.g., molecular systems if a generalized scaling parameter w = h~ s 
is introduced as a new dynamical variable [94]. Quantization yields bound 
states or resonances, Wk, for the scaling parameter. In scaling systems the 
semiclassical response function g sc (w) can be Fourier transformed easily to 
obtain the semiclassical trace of the propagator 

+oo 

C sc (s) = — f g s \w)e- lsw alw = £.A po 5 (s - s po ) . (3.32) 

271 -L po 
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The signal C sc (s) has <5-peaks at the positions of the classical periods (scaled 
actions) s = s po of periodic orbits and with peak heights (recurrence strengths) 
A po , i.e., C sc (s) is Gutzwiller's periodic orbit recurrence function. Consider 
now the quantum mechanical counterparts of g sc (w) and C sc (w) taken as the 
sums over the poles Wk of the Green's function, 

g ^( w ) = V d ± , (3.33) 



+00 



C* m (s) = — / g Cim (w)e- tsw dw = -iy j d k e- tWkS , (3.34) 



k 



with dk being the multiplicities of resonances, i.e., dk = 1 for non-degenerate 
states. In analogy with the calculation of Riemann zeros from Eq. (3.18) the 
frequencies, Wk, and amplitudes, dk, can be extracted by harmonic inversion 
of the signal C sc (s) after convoluting it with a Gaussian function, i.e., 

Cl c (s) = -L- ]T A P J S - S ^ 2 / 2 * 2 . (3.35) 

VZira po 



By adjusting C^(s) to the functional form of Eq. 3.34, the frequencies, w k , can 
be interpreted as the semiclassical approximation to the poles of the Green's 
function in (3.33). Note that the harmonic inversion method described in Ap- 
pendix A.l allows studying signals with complex frequencies Wk as well. For 
open systems the complex frequencies can be interpreted as semiclassical reso- 
nances. Note also that the Wk in general differ from the exact quantum eigen- 
values because Gutzwiller's trace formula (3.30) is an approximation, correct 
only to the lowest order in h. Therefore the diagonalization of small matrices 
(Eq. A. 7 in Appendix A.l) does not imply that the results of the periodic or- 
bit quantization are more "quantum" in any sense than those obtained, e.g., 
from a cycle expansion [9,58]. However, the harmonic inversion technique also 
allows the calculation of higher order % corrections to the periodic orbit sum, 
and we will return to this problem in Section 3.6. 



3.2.2 Semiclassical matrix elements 

The procedure described above can be generalized in a straightforward man- 
ner to the calculation of semiclassical diagonal matrix elements (ipk\A\ipk) of 
a smooth Hermitian operator A. In this case we start from the quantum me- 
chanical trace formula 

gTM = tr G+ A = £, ^ A ^\ (3.36) 

^ w - w k + i0 
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which has the same functional form as (3.33), but with d^ = (^|A|^) instead 
of dk = 1. For the quantum response function g ( ^ i {w) (Eq. 3.36) a semiclassical 
approximation has been derived in [128,129], which has the same form as 
Gutzwiller's trace formula (3.30) but with amplitudes 

A po = -i . Ap6 (3.37) 

v /|det(M po -/)| 



where M po is the monodromy matrix and yU po the Maslov index of the periodic 
orbit, and 



A p = J A(q(s),p(s))ds (3.38) 



is the classical integral of the observable A over one period s p of the primitive 
periodic orbit. Note that q(s) and p(s) are functions of the classical action 
instead of time for scaling systems [130]. Gutzwiller's trace formula for the 
density of states is obtained with A being the identity operator, i.e., A p = s p . 
When the semiclassical signal C sc (s) (Eq. 3.32) with amplitudes A po given 
by Eqs. 3.37 and 3.38 is analyzed with the method of harmonic inversion the 
frequencies and amplitudes obtained are the semiclassical approximations to 
the eigenvalues Wk and matrix elements dk = {ijJk\A\i(jk) , respectively. 



3.3 The three disk scattering system 



Let us consider the three disk scattering system which has served as a model 
for periodic orbit quantization by cycle expansion techniques [9-15]. For a 
brief introduction of the system, its symbolic dynamics and symmetries we 
refer the reader to Section 2.6.2. The starting point for the cycle expansion 
technique is to rewrite Gutzwiller's trace formula as a dynamical zeta function 
[131] 

v^na-^); ( 3 - 39 ) 



with p indicating the primitive periodic orbits, 

t p = -^= e i[SpW -^^ ] , (3.40) 
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Up the Maslov index, and A p the largest eigenvalue of the monodromy matrix. 
[The dynamical zeta function, Eqs. 3.39 and 3.40, is based on an approxi- 
mation. By contrast, the Gutzwiller-Voros zeta function [132,133] is exactly 
equivalent to Gutzwiller's trace formula. See Refs. [14,15] for more details.] 
The semiclassical eigenvalues w are obtained as zeros of the dynamical zeta 
function (3.39). However, the convergence problems of the infinite product 
(3.39) are similar to the convergence problems of Gutzwiller's periodic orbit 
sum. The basic observation for the three disk system with sufficiently large 
distance d was that the periodic orbit quantities of long orbits can be approx- 
imated by the periodic orbit quantities of short orbits, i.e., for two orbits with 
symbolic codes p\ and P2 we have 



l Jj PlP2 = A*pi ' A*P2 1 

and thus 

tpiP2 ~ %>1^P2 ■ V^-^^/ 



This implies that in the cycle expansion of the dynamical zeta function (3.39) 
the contributions of long periodic orbits are shadowed by the contributions of 
short orbits. E.g., for the three disk scattering system the cycle expansion up 
to cycle length n = 3 reads 

l/£ — 1 — to — ti — [toi — toti] — [tooi — to^oi + ton — toiti] — . . . (3.43) 



The terms to and t± are the fundamental contributions, while the terms in 
brackets are the curvature corrections, ordered by cycle length, and can rapidly 
decrease with increasing cycle length. With the cycle expansion (3.43) the dy- 
namical zeta function (3.39) is analytically continued to the physically impor- 
tant region below the real axis, where the resonances are located. For d = 6 
semiclassical resonances were calculated by application of the cycle expansion 
technique including all (symmetry reduced) periodic orbits up to cycle length 
n=13 [11,14,15]. 

We now demonstrate the usefulness of the harmonic inversion technique for 
the semiclassical quantization of the three disk scattering system. In contrast 
to the cycle expansion we will not make use of the symbolic code of the 
periodic orbits and the approximate relations (3.42) between the periodic orbit 
quantities. As in the previously discussed examples for billiard systems the 
scaled action s is given by the length L of orbits (s = L) and the quantized 
parameter is the absolute value of the wave vector k — |k| = \ / 2mE/h. For 
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the three disk system the periodic orbit signal C{L) reads 

_ °° T p-^^fJ-p 

C(L) = -i £ £ M p > S(L - rL p ) , (3.44) 

v r =\ J\det(Mr-I)\ 



with L p , M p , and \x v the geometrical length, monodromy matrix, and Maslov 
index of the symmetry reduced primitive periodic orbit p, respectively. The 
weight factors M p depend on the irreducible subspace (A ± , A 2 , and E), where 
resonances are calculated, and the symmetries of the orbits. Periodic orbits 
without symmetries, with rotational symmetry, and with symmetry under re- 
flection are characterized by e, C3, Cf , and a v , respectively. The weight factors 
M p are the same as for the cycle expansion, and are given in Table 7. For details 
of the symmetry decomposition see Ref. [134]. In the following we calculate 
resonances in the irreducible subspace Ai, i.e., M p = 1 in (3.44) for all orbits. 
We first apply harmonic inversion to the case R : d — 1 : 6 studied before. Fig. 
28a shows the periodic orbit recurrence function, i.e., the trace of the semiclas- 
sical propagator C SC (L). The groups with oscillating sign belong to periodic 
orbits with adjacent cycle lengths. To obtain a smooth function on an equidis- 
tant grid, which is required for the harmonic inversion method, the (^-functions 
in (3.44) have been convoluted with a Gaussian function of width a = 0.0015. 
As explained in Section 3.1.2 this does not change the underlying spectrum. 
The results of the harmonic inversion analysis of this signal are presented in 
Fig. 28b. The crosses in Fig. 28b represent semiclassical poles, for which the 
amplitudes dk are very close to 1, mostly within one percent. Because the am- 
plitudes converge much slower than the frequencies these resonance positions 
can be assumed to be very accurate within the semiclassical approximation. 
For some broad resonances marked by diamonds in Fig. 28b the dk deviate 
strongly from 1, within 5 to maximal 50 percent. It is not clear whether these 
strong deviations are due to numerical effects, such as convergence problems 
caused by too short a signal, or if they are a consequence of the semiclassi- 
cal approximation. A direct comparison between the semiclassical resonances 
k hl obtained by harmonic inversion of the periodic orbit sum, and results of 
Wirzba [14,15,113], i.e., the cycle expansion resonances k ce calculated up to 
12 th order in the curvature expansion of the Gutzwiller-Voros zeta function, 
and the exact quantum resonances k qm is given in Table 8 for resonances with 
Re k < 30 and in Table 9 for 120 < Re k < 132. Apart from a few resonances 
with large imaginary parts the differences \k hl — k cc \ between the semiclassical 
resonances obtained by harmonic inversion and cycle expansion are by several 
orders of magnitude smaller than the semiclassical error, i.e., the differences 

To compare the efficiency of the harmonic inversion and the cycle expansion 
method we calculated the A\ resonances of the three disk system with distance 
d = 6 from a short signal C(L) with L < 24, which includes the recurrences 
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of all periodic orbits with cycle length n < 5. The resonances obtained by 
harmonic inversion of the short signal are presented as squares in Fig. 29. The 
resonances of the two bands closest to the real axis qualitatively agree with the 
correct semiclassical resonances (crosses in Fig. 29). The accuracy is similar to 
the accuracy obtained by the cycle expansion up to 3 rd order in the curvature 
expansion [15], which includes all periodic orbits with cycle length n < 4. The 
reason for the somewhat higher efficiency of the cycle expansion compared to 
harmonic inversion is probably that the basic requirement (3.42) for the cycle 
expansion is a very good approximation at the large distance d = 6 between 
the disks. 

We now study the three disk scattering system with a short distance ratio 
d/R = 2.5. The signal C(L) is constructed from 356 periodic orbits with ge- 
ometrical length L < 7.5 (see Fig. 30a). For large L groups of orbits with 
the same cycle length of the symbolic code strongly overlap and cannot be 
recognized in Fig. 30a. Note that the signal contains complete sets of orbits 
up to topological length (cycle length) n = 9 only. The resonances obtained 
by harmonic inversion of the signal C(L) are presented in Fig. 30b. A compar- 
ison between the semiclassical resonances k hl obtained by harmonic inversion, 
the cycle expansion resonances k cc calculated up to 9 th order in the curva- 
ture expansion of the Gutzwiller-Voros zeta function, and the exact quantum 
resonances /c qm is given in Table 10. The cycle expansion and exact quan- 
tum calculations have been performed by Wirzba [113]. The results of both 
semiclassical methods are in good agreement, although a detailed comparison 
reveals that for some resonances towards the end of Table 10 the values of the 
cycle expansion are somewhat closer to the exact quantum mechanical results 
than those values obtained by harmonic inversion. 



3.4 Systems with mixed regular- chaotic dynamics 



In Section 3.1 we have applied harmonic inversion for the calculation of ze- 
ros of Riemann's zeta function as a mathematical model of a strongly chaotic 
bound system, and in Section 3.3 we have used harmonic inversion for the 
periodic orbit quantization of the three disk scattering system. Both systems 
have been solved with other especially designed methods, i.e., the Riemann 
zeta function with the help of the Riemann-Siegel formula and the three disk 
system by application of cycle expansion techniques. However, none of the 
special methods, which were designed to overcome the convergence problems 
of the semiclassical trace formula (see, e.g., Refs. [9,61,63,65]), has succeeded 
so far in correctly describing generic dynamical systems with mixed regular- 
chaotic phase spaces. In this section we want to demonstrate the universality 
of periodic orbit quantization by harmonic inversion by investigating generic 
dynamical systems. It will be the objective to contribute to solving the long- 
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standing problem of semiclassical quantization of nonintegrable systems in the 
mixed regular-chaotic regime. 

A first step towards the periodic orbit quantization of mixed systems has 
been done by Wintgen [78] on the hydrogen atom in a magnetic field. From 
Gutzwiller's truncated periodic orbit sum he obtained the smoothed part of the 
density of states presented in Fig. 31 at scaled energy E = —0.2 which is in the 
mixed regular-chaotic regime. The resolution of the smoothed spectra depends 
on the cutoff value of the periodic orbit sum which is s max /2vr = 1.33 and 
s m ax/27T = 3.00 in Fig. 31a and 31b, respectively. Nine approximate eigenvalues 
have been estimated from the low-resolution truncated periodic orbit sum in 
Fig. 31b (for comparison the exact quantum eigenvalues are marked by vertical 
bars in Fig. 31), but neither the accuracy nor the number of semiclassical 
eigenvalues could be improved by increasing the number of periodic orbits 
because of the non-convergence property of Gutzwiller's trace formula. The 
results of Ref. [78] therefore clearly demonstrate that methods to overcome 
the convergence problems of the semiclassical periodic orbit sum are of crucial 
importance to obtain highly resolved semiclassical spectra in the mixed regime. 
Note also that Gutzwiller's trace formula is not valid for non-isolated periodic 
orbits on invariant tori in the regular part of the phase space. This fact has 
not been considered in Ref. [78]. 

The limiting cases of mixed systems are strongly chaotic and integrable sys- 
tems. As has been proven by Gutzwiller [8], for systems with complete chaotic 
(hyperbolic) classical dynamics the density of states can be expressed as an 
infinite sum over all (isolated) periodic orbits. On the other extreme of com- 
plete integrability, it is well known that the semiclassical energy values can 
be obtained by EBK torus quantization [1] . This requires the knowledge of all 
the constants of motion, which are not normally given in explicit form, and 
therefore practical EBK quantization based on the direct or indirect numerical 
construction of the constants of motion turns out to be a formidable task [135]. 
As an alternative, EBK quantization was recast as a sum over all periodic or- 
bits of a given topology on respective tori by Berry and Tabor [25,26]. The 
Berry- Tabor formula circumvents the numerical construction of the constants 
of motion but usually suffers from the convergence problems of the infinite 
periodic orbit sum. 

The extension of the Berry- Tabor formula into the near-integrable (KAM) 
regime was outlined by Ozorio de Almeida [41] and elaborated, at different 
levels of refinement, by Tomsovic et al. [27] and Ullmo et al. [28]. These authors 
noted that in the near-integrable regime, according to the Poincare-Birkhoff 
theorem, two periodic orbits survive the destruction of a rational torus with 
similar actions, one stable and one hyperbolic unstable, and worked out the 
ensuing modifications of the Berry- Tabor formula. In this section we go one 
step further by noting that, with increasing perturbation, the stable orbit 
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turns into an inverse hyperbolic one representing, together with its unstable 
companion with similar action, a remnant torus. We include the contributions 
of these pairs of inverse hyperbolic and hyperbolic orbits in the Berry- Tabor 
formula and demonstrate for a system with mixed regular-chaotic dynamics 
that this procedure yields excellent results even in the deep mixed regular- 
chaotic regime. As in Ref. [78] we choose the hydrogen atom in a magnetic 
field, which is a real physical system and has served extensively as a prototype 
for the investigation of "quantum chaos" [29-31]. 

The fundamental obstacle bedeviling the semiclassical quantization of systems 
with mixed regular-chaotic dynamics is that the periodic orbits are neither 
sufficiently isolated, as is required for Gutzwiller's trace formula [8], nor are 
they part of invariant tori, as is necessary for the Berry- Tabor formula [25,26]. 
However, as will become clear below, it is the Berry- Tabor formula which lends 
itself in a natural way for an extension of periodic orbit quantization to mixed 
systems. As previously, we consider scaling systems where the classical action 
S scales as S = sw with w = h~$ the scaling parameter and s the scaled 
action. For scaling systems with two degrees of freedom, which we will focus 
on, the Berry- Tabor formula for the fluctuating part of the level density reads 

e ( w ) = I Re V W '! 2sm e <(«M W -i™-5) ( 3.4 5) 

* mM 2 3/2 K|V2 



with M = (Mi, M 2 ) pairs of integers specifying the individual periodic orbits 
on the tori (numbers of rotations per period, M 2 /Ml rational), and «m and ??m 
the scaled action and Maslov index of the periodic orbit M. The function g% 
in (3.45) is obtained by inverting the Hamiltonian, expressed in terms of the 
actions (ii,/ 2 ) of the corresponding torus, with respect to I2, viz. H(Ii,l2 = 
9e{I\J) — E [136]. The calculation of g" E from the actions (ii, I2) can be rather 
laborious even for integrable and near-integrable systems, and, by definition, 
becomes impossible for mixed systems in the chaotic part of the phase space. 
Here we will adopt the method of Refs. [27,28] and calculate g E , for given 
M = (yUi,/i 2 ), with (/!i,/i 2 ) coprime integers specifying the primitive periodic 
orbit, directly from the parameters of the two periodic orbits (stable (s) and 
hyperbolic unstable (h)) that survive the destruction of the rational torus M, 
viz. 

9e = —T^T I / + = I , (3.46) 

wfiiAs Wdet(M s - /) J- det (M h - /) 



with 



As = ^(s h -s s ) (3.47) 
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the difference of the scaled actions, and M s and M^ the monodromy matrices 
of the two orbits. The action sm in (3.45) is to replaced with the mean action 

s = -{s h + s s ). (3.48) 



Eq. 3.46 is an approximation which becomes exact in the limit of an integrable 
system. 

It is a characteristic feature of systems with mixed regular-chaotic dynamics 
that with increasing nonintegrability the stable orbits turn into inverse hy- 
perbolic unstable orbits in the chaotic part of the phase space. These orbits, 
although embedded in the fully chaotic part of phase space, are remnants of 
broken tori. It is therefore natural to assume that Eqs. 3.45 and 3.46 can even 
be applied when these pairs of inverse hyperbolic and hyperbolic orbits are 
taken into account, i.e., more deeply in the mixed regular-chaotic regime. 

It should be noted that the difference As between the actions of the two orbits 
is normally still small, and it is therefore more appropriate to start from 
the Berry- Tabor formula for semiclassical quantization in that regime than 
from Gutzwiller's trace formula, which assumes well-isolated periodic orbits. 
It is also important to note that the Berry- Tabor formula does not require 
an extensive numerical periodic orbit search. The periodic orbit parameters 
s/M 2 and g" E are smooth functions of the rotation number M 2 /Mi, and can 
be obtained for arbitrary periodic orbits with coprime integers (Mi,M2) by 
interpolation between "simple" rational numbers M 2 jM\. 



3.4-1 Hydrogen atom in a magnetic field 

We now demonstrate the high quality of the extension of Eqs. 3.45 and 3.46 
to pairs of inverse hyperbolic and hyperbolic periodic orbits for a physical sys- 
tem that undergoes a transition from regularity to chaos, namely the hydrogen 
atom in a magnetic field. This is a scaling system, with w = 7 -1 / 3 = h~g the 
scaling parameter and 7 = -8/(2.35 x 10 5 T) the magnetic field strength in 
atomic units. The classical Hamiltonian and the scaling procedure have been 
discussed in Section 2.1. At low scaled energies E = E^~ 2 ^ < —0.6 a Poincare 
surface of section analysis [30] of the classical Hamiltonian (Eq. 2.11) exhibits 
two different torus structures related to a "rotator" and "vibrator" type mo- 
tion. The separatrix between these tori is destroyed at a scaled energy of 
E ?» —0.6, and the chaotic region around the separatrix grows with increasing 
energy. At E — —0.127 the classical phase space becomes completely chaotic. 
We investigate the system at scaled energy E = —0.4, where about 40% of 
the classical phase space volume is chaotic (see inset in Fig. 32), i.e. well in 
the region of mixed dynamics. We use 8 pairs of periodic orbits to describe 
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the rotator type motion in both the regular and chaotic region. The results 
for the periodic orbit parameters s/27rM 2 and g E are presented as solid lines 
in Fig. 32. The squares on the solid lines mark parameters obtained by pairs 
of stable and unstable periodic orbits in the regular region of the phase space. 
The diamonds mark parameters obtained by pairs of two unstable (inverse hy- 
perbolic and hyperbolic) periodic orbits in the chaotic region of phase space. 
The cutoff is related to the winding angle = 1.278 of the fixed point of the 
rotator type motion, i.e., the orbit perpendicular to the magnetic field axis, 
(M 2 /M 1 ) cutoff = 7r/0 = 2.458. The solid lines have been obtained by spline 
interpolation of the data points. In the same way the periodic orbit parame- 
ters for the vibrator type motion have been obtained from 11 pairs of periodic 
orbits (see the dashed lines in Fig. 32). The cutoff at M 2 /Mi = 7r/0 = 1.158 
is related to the winding angle = 2.714 of the fixed point of the vibrator 
type motion, i.e., the orbit parallel to the field axis. 

With the data of Fig. 32 we have all the ingredients at hand to calculate 
the semiclassical density of states g(w) in Eq. 3.45. The periodic orbit sum in- 
cludes for both the rotator and vibrator type motion the orbits with M 2 /Mi > 
(M 2 /M 1 ) cuto g. For each orbit the action and the function g" E is obtained from 
the spline interpolations. The Maslov indices are rj-^i = 4M 2 — Mi for the 
rotator and t/m — 4M 2 + 1M\ — 1 for the vibrator type orbits. However, the 
problem is to extract the semiclassical eigenenergies from Eq. 3.45 because the 
periodic orbit sum does not converge. We now adopt the method discussed 
in Section 3.2 and adjust the semiclassical recurrence signal, i.e., the Fourier 
transform of the weighted density of states w~ 1 ^ 2 g(w) (Eq. 3.45) 

C sc (s) = J2A M 5(s-s M ), (3.49) 

M 



with the amplitudes being determined exclusively by periodic orbit quantities, 

Am = J M - e-'?™ , (3.50) 

M 2 3/ Vi;l 1/2 



to the functional form of its quantum mechanical analogue 

C* a (s) = -iJ2dke- iWka , (3.51) 

k 



where the Wk are the quantum eigenvalues of the scaling parameter, and the dk 
are the multiplicities of the eigenvalues (dk = 1 for nondegenerate states). The 
frequencies obtained from this procedure are interpreted as the semiclassical 
eigenvalues w^. The technique used to adjust Eq. (3.49) to the functional form 
of Eq. (3.51) is harmonic inversion [82,83]. 
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For the hydrogen atom in a magnetic field part of the semiclassical recurrence 
signal C sc (s) at scaled energy E = —0.4 is presented in Fig. 33. The solid and 
dashed peaks mark the recurrencies of the rotator and vibrator type orbits, 
respectively. Note that C sc (s) can be easily calculated even for long periods s 
with the help of the spline interpolation functions in Fig. 32. By contrast, the 
construction of the recurrence signal for Gutzwiller's trace formula usually 
requires an exponentially increasing effort for the numerical periodic orbit 
search with growing period. 

We have analyzed C sc (s) by the harmonic inversion technique in the region 
< s/2n < 200. The resulting semiclassical spectrum of the lowest 106 states 
with eigenvalues w < 20 is shown in the upper part of Fig. 34a. For graphical 
purposes the spectrum is presented as a function of the squared scaling param- 
eter w 2 , which is equivalent to unfolding the spectrum to constant mean level 
spacing. For comparison the lower part of Fig. 34a shows the exact quantum 
spectrum. The semiclassical and quantum spectrum are seen to be in excellent 
agreement, and deviations are less than the stick widths for nearly all states. 
The distribution P(d) of the semiclassical error with d = (u> qm — w sc )/Aw av 
the error in units of the mean level spacing, Aw av = 1.937/w, is presented in 
Fig. 34b. For most levels the semiclassical error is less than 4% of the mean 
level spacing, which is typical for a system with two degrees of freedom [110]. 

The accuracy of the results presented in Fig. 34 seems to be surprising for 
two reasons. First, we have not exploited the mean staircase function J\f(w), 
i.e., the number of eigenvalues Wk with Wk < w, which is a basic require- 
ment of some other semiclassical quantization techniques for bound chaotic 
systems [65,61]. Second, as mentioned before, Eq. 3.46 has been derived for 
near-integrable systems, and is only an approximation, in particular, for mixed 
systems. We have not taken into account any more refined extensions of the 
Berry- Tabor formula (3.45) as discussed, e.g., in Refs. [27,28]. The answer to 
the second point is that the splitting of scaled actions of the periodic orbit 
pairs used in Fig. 32 does not exceed As = 0.022, and therefore for states 
with w < 20 the phase shift between the two periodic orbit contributions is 
wAs = 0.44, at most. For small phase shifts the extension of the Berry- Tabor 
formula to near-integrable systems results in a damping of the amplitudes of 
the periodic orbit recurrence signal in Fig. 33 but seems not to effect the 
frequencies, i.e., the semiclassical eigenvalues Wk obtained by the harmonic 
inversion of the function C sc (s). 

To summarize, in this section we have presented a solution to the fundamental 
problem of semiclassical quantization of nonintegrable systems in the mixed 
regular-chaotic regime. We have demonstrated the excellent quality of our 
procedure for the hydrogen atom in a magnetic field at a scaled energy E = 
—0.4, where about 40% of the phase space volume is chaotic. The lowest 106 
semiclassical and quantum eigenenergies have been shown to agree within a 
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few percent of the mean level spacings. The same method can be applied 
straightforwardly to other systems with mixed dynamics. 



3.5 Harmonic inversion of cross- correlated periodic orbit sums 



In the previous sections we have introduced harmonic inversion of semiclassi- 
cal signals as a powerful and universal technique for the problem of periodic 
orbit quantization in that it does not depend on special properties of the 
system such as being bound and ergodic, or the existence of a symbolic dy- 
namics. The method only requires the knowledge of periodic orbits and their 
physical quantities up to a certain maximum period, which depends on the 
average local density of states. Unfortunately, this method is not free of the 
general drawback of most semiclassical approaches, which suffer from a rapid 
proliferation of periodic orbits with their period, which in turn requires an 
enormous number of orbits to be taken into account. As discussed in Section 
3.1.4 the required signal length for harmonic inversion is at least two times 
the Heisenberg period. This is by a factor of four longer than the signal length 
required for the application of the Riemann-Siegel look-alike formula [59-61]. 
In chaotic systems, where the number of periodic orbits grows exponentially 
with the period, the reduction of the required signal length by a factor of four 
implies that, e.g., instead of one million periodic orbits a reduced set of about 
32 orbits is sufficient for the semiclassical quantization. This example clearly 
illustrates that a shortening of the required signal length is highly desirable. 

In this section we want to introduce harmonic inversion of cross-correlated 
periodic orbit sums as a method to reduce the required amount of periodic 
orbit data [85] . The idea is that the informational content of a D x D cross- 
correlated time signal is increased roughly by a factor of D as compared to 
a 1 x 1 signal. The cross-correlated signal is constructed by introducing a set 
of D smooth and linearly independent operators. Numerically, the harmonic 
inversion of the cross-correlated periodic orbit sum is based on an extension of 
the filter-diagonalization method to the case of time cross-correlation functions 
[68,137,138]. This extended method provides highly resolved spectra even in 
situations of nearly degenerate states, as well as the diagonal matrix elements 
for the set of operators chosen. The power of the method will be demonstrated 
in Section 3.7 for the circle billiard, as an example of a completely integrable 
system. 

Consider a quantum Hamiltonian H whose eigenvalues are w n and eigenstates 
\n). [As previously, we consider scaling systems with w the scaling parameter, 
which is not necessarily the energy, and with linear scaling of the classical 
action, S po = ws po .} We introduce a cross-correlated response function (a, (3 = 
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1,2,. ..,D) 



^— v "an"/3n 

^^ W — W n + iO 



where b an and b@ n are the diagonal matrix elements of two operators A a and 
A/3, respectively, i.e. 

b an = (n\A a \n) . (3.53) 



Later we have to find a semiclassical approximation to Eq. 3.52. In this context 
it is important to note that g a p(w) can only be written as a trace formula, 



g. 

with the Green Function 

G + (w) 



a/3 (w) =tr {A a G + (w)Ap} 



w - H + i0 ' 



if either A a or Ap commutes with H . The weighted density of states is given 
by 

Q a p{w) = lrag af3 (w) . (3.54) 

7T 



Let us assume that the semiclassical approximation g^piw) to the quantum 
expression (3.52) is given. The general procedure of harmonic inversion as 
described in Section 3.2 would then be to adjust the Fourier transform of 
g s £Jw) to the functional form of the quantum expression 

+oo 

°^ {S) = 27 / 3T^y~ lsw dw = -i £ b an b Pn e-^ s . (3.55) 



The conventional harmonic inversion problem is formulated as a nonlinear fit 
of the signal C(s) by the sum of sinusoidal terms (see Section 2.2 and Appendix 

A.l), 

C(s) = J2d n e- tw " s , 
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with the set of, in general, complex variational parameters {w n , d n = —ib an b/3 n }. 
As already discussed in Section 3.1.4, simple information theoretical consid- 
erations then yield an estimate for the required signal length, s max ~ 4iiq(w), 
for poles w n < w which can be extracted by this method. When a periodic 
orbit approximation of the quantum signal C(s) is used, this estimate results 
sometimes in a very unfavorable scaling because of a rapid (exponential for 
chaotic systems) proliferation of periodic orbits with increasing period. Let us 
consider now a generalized harmonic inversion problem, which assumes that 
the whole s-dependent D x D signal C a p(s) is adjusted simultaneously to the 
form of Eq. 3.55, with b an and w n being the variational parameters. The ad- 
vantage of using the cross-correlation approach [68,137,138] is based on the 
simple argument that the total amount of independent information contained 
in the D x D signal is D(D + 1) multiplied by the length of the signal, while 
the total number of unknowns (here b an and w n ) is (D + 1) times the total 
number of poles w n . Therefore the informational content of the D x D sig- 
nal per unknown parameter is increased (compared to the case of Eq. 3.32) 
by a factor of D. [Of course, this scaling holds only approximately and for 
sufficiently small numbers D of operators A a chosen.] Thus we have the re- 
sult that, provided we are able to obtain a periodic orbit approximation for 
C a /3{s), with this procedure we can extract more information from the same 
set of periodic orbits. 

We now have to find a semiclassical approximation for C a p(s) in Eq. 3.55. The 
problem has been solved in the literature only for special cases, i.e., if one or 
both operators A a and Ap are the identity or, somewhat more general, if at 
least one of the operators commutes with the Hamiltonian H . For the identity 
operator A± — I the element C^_(s) is the Fourier transform of Gutzwiller's 
trace formula [7,8] for chaotic systems, and of the Berry- Tabor formula [25,26] 
for regular systems, i.e. 

CTi(s) = J2AoHs-s po ) , (3.56) 

po 



where s po are the periods of the orbits and A po the amplitudes (recurrence 
strengths) of the periodic orbit contributions including phase information. For 
A\ — I and an arbitrary smooth operator A a the elements C^^s) are obtained 
from a semiclassical approximation to the generalized trace formula 

g al (w)=ti{G + (w)A a } , (3.57) 



which has been investigated in detail in Refs. [128,129]. The result is that the 
amplitudes A po in (3.56) have to be multiplied by the classical average of the 
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observable A a along the periodic orbit, 

Spo 

a a , po = — / A a (q(s),p(s))ds , (3.5* 



po 



with A a (q, p) the Wigner transform of the operator A a , i.e., the signal C^i(s) 
reads 

Cl c i{s) = J2 a*,poA po S (s - Spo) • (3.59) 

po 



If at least one of the operators A a and Ap commutes with H, Eq. 3.52 can 
still be written as a trace formula, 

g a0 (w) = tr {A a G + (w)Ap} = tr {G + (w) (A a Ap)} , 



and Eq. 3.59 can be applied to the product A a Ap. However, we do not want 
to restrict the operators to those commuting with H, which obviously would 
be a severe restriction especially for chaotic systems, and the problem is now 
to find a semiclassical approximation to Eq. 3.55 for the general case of two 
arbitrary smooth operators A a and Ap. A reasonable assumption is that the 
amplitudes A po in (3.56) have to be multiplied by the product of the classical 
averages, a a ^ ap^ 0l of these two observables, i.e. 

CqM 5 ) = H a a,V° a P,V°A po 5 (S - Spo) • (3.60) 

po 

Although no rigorous mathematical proof of Eq. 3.60 will be given here, we 
have strong numerical evidence, from the high resolution analysis of quantum 
spectra, part of which will be given below, that the conjecture of Eq. 3.60 
is correct. Details on semiclassical non-trace type formulae like Eq. 3.60 are 
given in Ref. [139]. Eq. 3.60 is the starting point for the following application 
of harmonic inversion of cross-correlation functions. Note that all quantities 
in (3.60) are obtained from the classical periodic orbits. 

The idea of periodic orbit quantization by harmonic inversion of cross-correlated 
periodic orbit sums is to fit the semiclassical functions C^pis) given in a finite 
range < s < s max to the functional form of the quantum expression (3.55). 
As for the harmonic inversion of a one-dimensional signal (see Section 3.2) 
the frequencies of the harmonic inversion analysis are then identified with the 
semiclassical eigenvalues w n . The amplitudes b an are identified with the semi- 
classical approximations to the diagonal matrix elements (n\A a \n). Here we 
only give a brief description how the harmonic inversion method is extended 
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to cross-correlation functions. The details of the numerical procedure of solv- 
ing the generalized harmonic inversion problem (3.55) are presented in Refs. 
[68,137,138] and in Appendix A. 2. As for the harmonic inversion of a single 
function the idea is to recast the nonlinear fit problem as a linear algebraic 
problem [68]. This is done by associating the signal C a p(s) (to be inverted) 
with a time cross-correlation function between an initial state $ a and a final 
state $g, 

C a p(s) = (® p ,e- iCls <f> a ) , (3.61) 



where the fictitious quantum dynamical system is described by an effective 
Hamiltonian Cl. The latter is defined implicitly by relating its spectrum to the 
set of unknown spectral parameters w n and b an . Diagonalization of Cl would 
yield the desired w n and b an . This is done by introducing an appropriate basis 
set in which the matrix elements of fl are available only in terms of the known 
signals C a p(s). The Hamiltonian Cl is assumed to be complex symmetric even 
in the case of a bound system. This makes the harmonic inversion stable with 
respect to "noise" due to the imperfections of the semiclassical approximation. 
The most efficient numerical and practical implementation of the harmonic 
inversion method with all relevant formulae can be found in Refs. [137,138] 
and Appendix A. 2. 

The method of harmonic inversion of cross-correlated periodic orbit sums will 
be applied in Section 3.7 to the circle billiard. As will be shown, for a given 
number of periodic orbits the accuracy of semiclassical spectra can be signif- 
icantly improved with the help of the cross-correlation approach, or, alterna- 
tively, spectra with similar accuracy can be obtained from a periodic orbit 
cross-correlation signal with significantly reduced signal length. 



3.6 h expansion for the periodic orbit quantization by harmonic inversion 



Semiclassical spectra can be obtained for both regular and chaotic systems in 
terms of the periodic orbits of the system. For chaotic dynamics the semiclas- 
sical trace formula was derived by Gutzwiller [7,8], and for integrable systems 
the Berry- Tabor formula [25,26] is well known to be precisely equivalent to 
the EBK torus quantization [1-3]. However, as already has been discussed 
in Section 2.6, the semiclassical trace formulae are exact only in exceptional 
cases, e.g., the geodesic motion on the constant negative curvature surface. In 
general, they are just the leading order terms of an infinite series in powers 
of the Planck constant and the accuracy of semiclassical quantization is still 
an object of intense investigation [109,110,140]. Methods for the calculation of 
the higher order periodic orbit contributions were developed in Refs. [22-24]. 
In Section 2.6 we have demonstrated how the periodic orbit quantities of the 
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h expanded trace formula can be extracted from the quantum and semiclas- 
sical spectra. It is an even more fundamental problem to obtain semiclassical 
eigenenergies beyond the Gutzwiller and Berry- Tabor approximation directly 
from the H expanded periodic orbit sum. Note that the H expansion of the 
periodic orbit sum does not solve the general problem of the construction of 
the analytic continuation of the trace formula, which is already a fundamental 
problem when only the leading order terms of the h expansion is considered. 
Up to now the h expansion for periodic orbit quantization is restricted to sys- 
tems with known symbolic dynamics, like the three disk scattering problem, 
where cycle expansion techniques can be applied [23,24], and semiclassical 
eigenenergies beyond the Gutzwiller and Berry- Tabor approximation cannot 
be calculated, e.g., for bound systems with the help of Riemann-Siegel type 
formulae [60,61] or surface of section techniques [62,63]. In this section we ex- 
tend the method of periodic orbit quantization by harmonic inversion to the 
analysis of the % expansion of the periodic orbit sum. The accuracy of semi- 
classical eigenvalues can be improved by one to several orders of magnitude, 
as will be shown in Section 3.7 by way of example of the circle billiard. 

As in Section 2.6 we consider systems with a scaling property, i.e., where the 
classical action scales as S po = ws po , and the scaling parameter w = h^ 
plays the role of an inverse effective Planck constant. The h expansion of the 
periodic orbit sum is given (see Eq. 2.65) as a power series in w' 1 , 

OO OO 1 

9H = E 9n(w) = E — E A$e is ^ . (3.62) 

ra=0 n=0 W po 



The complex amplitudes A p n J of the n th order periodic orbit contributions 
include the phase information from the Maslov indices. For periodic orbit 
quantization the zeroth order contributions A p ^ are usually considered only. 
The Fourier transform of the principal periodic orbit sum 

+oo 

C (s) = — I g (w)e^ sw dw = E4°o 5 ( s " s p°) ( 3 - 63 ) 

27r io 



is adjusted by application of the harmonic inversion technique (see Section 
3.2) to the functional form of the exact quantum expression 

1 + r° H 

C(s) = — / Y e~ lws dw = -iY d k e- lw * s , (3.64) 

1 ; 2tt J ^ w - w k + tO Y 



with {wk, dk} the eigenvalues and multiplicities. The frequencies Wkfl obtained 
by harmonic inversion of Eq. 3.63 are the zeroth order h approximation to 
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the semiclassical eigenvalues. We will now demonstrate how the higher or- 
der correction terms to the semiclassical eigenvalues can be extracted from 
the periodic orbit sum (3.62). We first remark that the asymptotic expansion 
(3.62) of the semiclassical response function suffers, for n > 1, from the sin- 
gularities at w — 0, and it is therefore not appropriate to harmonically invert 
the Fourier transform of (3.62), although the Fourier transform formally exists. 
This means that the method of periodic orbit quantization by harmonic inver- 
sion cannot straightforwardly be extended to the h expansion of the periodic 
orbit sum. Instead we will calculate the correction terms to the semiclassical 
eigenvalues separately, order by order, as described in the following. 

Let us assume that the (n — l) st order approximations w k , n -i to the semi- 
classical eigenvalues are already obtained and the w k>n are to be calculated. 
The difference between the two subsequent approximations to the quantum 
mechanical response function reads [86] 



k \ 



dh di 



W - W k „ + iO W - Wkn-l + i® , 



^ d k Aw kyn 
k {w - w Kn + iO) 2 

with w k>n = (w kjn + w fc ,n-i)/2 and Aw k>n = w k>n - w k , n -i- Integration of Eq. 
3.65 and multiplication by w n yields 

Q n {w) = w / g n (w)dw = > — - , (3.66) 

J u w - w kn + iO 



k 



which has the functional form of a quantum mechanical response function but 
with residues proportional to the n th order corrections Aw ktU to the semiclas- 
sical eigenvalues. The semiclassical approximation to Eq. 3.66 is obtained from 
the term g n (w) in the periodic orbit sum (3.62) by integration and multipli- 
cation by w n , i.e. 

G n (w) = w n fg n (w)dw = -iJ2—A$e iws *° + o(-) . (3.67) 
We can now Fourier transform both Eqs. 3.66 and 3.67, and obtain (n > 1) 

C n {s) = — f g n (w)e^ ws dw = i £ d k (w k ) n Aw k , n e- iw " s (3.68) 



- DO 



= -*£ — A$S(s-* P o). (3.69) 

po S po 
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Eqs. 3.68 and 3.69 are the main result of this section. They imply that the 
h expansion of the semiclassical eigenvalues can be obtained, order by order, 
by harmonic inversion (h.i.) of the periodic orbit signal in Eq. 3.69 to the 
functional form of Eq. 3.68. The frequencies of the periodic orbit signal (3.69) 
are the semiclassical eigenvalues Wk- Note that the accuracy of the semiclassical 
eigenvalues does not necessarily increase with increasing order n. We indicate 
this in Eq. 3.68 by omitting the index n at the eigenvalues Wk- The corrections 
Awk t n to the eigenvalues are obtained from the amplitudes, dk(wk) n Awk t n, of 
the periodic orbit signal. 

The method requires as input the periodic orbits of the classical system up to 
a maximum period (scaled action), s max , determined by the average density 
of states [82,83]. The amplitudes A^ are obtained from Gutzwiller's trace 
formula [7,8] and the Berry- Tabor formula [25,26] for chaotic and regular sys- 
tems, respectively. For the next order correction A$ explicit formulae were 
derived by Gaspard and Alonso for chaotic systems with smooth potentials 
[22] and in Refs. [23,24] for billiards. With appropriate modifications [87,88] 
the formulae can be used for regular systems as well. As an example we inves- 
tigate the ft, expansion of the periodic orbit sum for the circle billiard in the 
next section. 



3.7 The circle billiard 



In Section 3.5 (see also [85]) we have introduced harmonic inversion of cross- 
correlated periodic orbit sums as a method to significantly reduce the required 
number of periodic orbits for semiclassical quantization, and in Section 3.6 (see 
also [86]) we have discussed the ft expansion of the periodic orbit sum and the 
calculation of semiclassical eigenenergies beyond the Gutzwiller [8] and Berry- 
Tabor [25,26] approximation. We now demonstrate both methods by way of 
example of the circle billiard. The circle billiard is a regular system and has 
been chosen here for the following reasons. 

(1) The nearest neighbor level statistics of integrable systems is a Poisson 
distribution, with a high probability for nearly degenerate states. The 
conventional method for periodic orbit quantization by harmonic inver- 
sion requires very long signals to resolve the nearly degenerate states. 
We will demonstrate the power of harmonic inversion of cross-correlated 
periodic orbit sums by fully resolving those nearly degenerate states with 
a significantly reduced set of orbits. 

(2) All relevant physical quantities, i.e., the quantum and semiclassical eigen- 
energies, the matrix elements of operators, the periodic orbits and their 
zeroth and first order amplitudes of the ft expanded periodic orbit sum, 
and the periodic orbit averages of classical observables can easily be ob- 



tained. 
(3) The semiclassical quantization of the circle billiard as an example of an 
integrable system demonstrates the universality and wide applicability 
of periodic orbit quantization by harmonic inversion, i.e., the method 
is not restricted to systems with hyperbolic dynamics like, e.g., pinball 
scattering. 

The circle billiard has already been introduced in Section 2.6.1. The exact 
quantum mechanical eigenvalues E = h 2 k 2 /2M are given as zeros of Bessel 
functions J\ m \(kR) = 0, where m is the angular momentum quantum number 
and R, the radius of the circle. In the following we choose R — 1. The semi- 
classical eigenvalues are obtained from the EBK quantization condition, Eq. 
2.71, 

\m\ / 3 



kRJl — (m/kR) 2 — \m\ arccos - — = tx [n -\ — 
v ' ' kR V 4 



and the h expansion of the periodic orbit sum reads (see Eq. 2.73) 



-I 1 OO OO 1 

V 1 V ff n=0 n=0 ^ po 



with (see Eqs. 2.74 and 2.75) 



, ( o)_ j![Ql_ p -»(f/ipo+5) 

po V 2 m y 

,m 1 / 5 — 2 sin 2 7 _, Y2 .. ^ 

4o = oV^ q . 3/2 e t{ ^°^ ) - 

I 3sm ' 7 

The angle 7 is defined as 7 = irm^/rrir, with m<$> = 1,2,... the number of turns 
of the periodic orbit around the origin, and m r = 2m ( / ) , 2m^ + l, . . . the number 
of reflections at the boundary of the circle. £ po = 2m r sin 7 and /i po = 3m r are 
the geometrical length and Maslov index of the orbits, respectively. 



3.7.1 Harmonic inversion of the cross-correlated periodic orbit sum 

We now calculate the semiclassical eigenenergies of the circle billiard by har- 
monic inversion of the cross-correlated periodic orbit sum (3.60) with A po = 
A p ^ (Eq. 2.74) the amplitudes of the Berry- Tabor formula [25,26], i.e., the low- 
est order H approximation. To construct the periodic orbit cross- correlation 
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signal C^g(£) we choose three different operators, 



the identity, 



A 2 = r 



the distance from the origin, and 



A 3 = (L/kf 



the square of the scaled angular momentum. For these operators the classical 
weights a^po (Eq. 3.58) are obtained as 



a l,po — 1 

1 / cos 7 
- 1 + 

2 y tan 7 

a3, P o = cos 2 7 . 



a 2, P o = - ( 1 + — arsinh tan 7 ) (3.70) 



Once all the ingredients of Eq. 3.60 for the circle billiard are available, the 3x3 
periodic orbit cross-correlation signal C^J£) can easily be constructed and in- 
verted by the generalized filter-diagonalization method. Results obtained from 
the periodic orbits with maximum length s max = 100 are presented in Fig. 35. 
Fig. 35a is part of the density of states, g(k), Figs. 35b and 35c are the den- 
sity of states weighted with the diagonal matrix elements of the operators 
A = r and A = L 2 , respectively. The squares are the results from the har- 
monic inversion of the periodic orbit cross-correlation signals. For comparison 
the crosses mark the matrix elements obtained by exact quantum calculations 
at positions k EBK obtained from the EBK quantization condition (2.71). In 
this section we do not compare with the exact zeros of the Bessel functions 
because Eq. 3.60 is correct only to first order in H and thus the harmonic 
inversion of C^Js) cannot provide the exact quantum mechanical eigenval- 
ues. The calculation of eigenenergies beyond the Berry- Tabor approximation 
will be discussed in Section 3.7.2. However, the perfect agreement between 
the eigenvalues k m obtained by harmonic inversion and the EBK eigenvalues 
k EBK is remarkable, and this is even true for nearly degenerate states marked 
by arrows in Fig. 35a. The eigenvalues of some nearly degenerate states are 
presented in Table 11. It is important to emphasize that these states with level 
splittings of, e.g., Ak = 6 x 10~ 4 cannot be resolved by the originally proposed 
method of periodic orbit quantization by harmonic inversion (see Section 3.2) 
with a periodic orbit signal length s max = 100. To resolve the two levels at 
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k m 11.049 (see Table 11) a signal length of at least s max ~ 500 is required 
if a single periodic orbit function C sc (s) is used instead of a cross-correlation 
function. The method presented in Section 3.5 can therefore be used to sig- 
nificantly reduce the required signal length and thus the required number of 
periodic orbits for periodic orbit quantization by harmonic inversion. As such 
the part of the spectrum shown in Fig. 35 can even be resolved, apart from the 
splittings of the nearly degenerate states marked by the arrows, from a short 
cross-correlation signal with s max = 30, which is about the Heisenberg period 
Sh = 27r^(A;), i.e. half of the signal length required for the harmonic inversion 
of a 1 x 1 signal. With five operators and a 5 x 5 cross-correlation signal highly 
excited states around k = 130 have even been obtained with a signal length 
Smax ~ 0.7sh [87,88], which is close to the signal length s max m 0.5sh required 
for the Riemann-Siegel type quantization [59-61]. The reduction of the signal 
length is especially important if the periodic orbit parameters are not given 
analytically, as in our example of the circle billiard, but must be obtained from 
a numerical periodic orbit search. How small can s max get as one uses more 
and more operators in the method? It might be that half of the Heisenberg 
period is a fundamental barrier for bound systems with chaotic dynamics in 
analogy to the Riemann-Siegel formula [61] while for regular systems an even 
further reduction of the signal length should in principle be possible. However, 
further investigations are necessary to clarify this point. 



3.7.2 Periodic orbit quantization beyond the Berry-Tabor approximation 

The semiclassical eigenvalues obtained by harmonic inversion of a cross-corre- 
lated or a sufficiently long single signal are in excellent agreement with the 
results of the EBK torus quantization, Eq. 2.71. However, they deviate from 
the exact quantum mechanical eigenenergies, i.e., the zeros of the Bessel func- 
tions because the Berry- Tabor formula [25,26] is only the lowest order Ti ap- 
proximation of the periodic orbit sum. We now demonstrate the % expansion 
of the periodic orbit sum and apply the technique discussed in Section 3.6 to 
the circle billiard. 

The first order corrections to the semiclassical eigenvalues, Ak = k^ — k^ 
are obtained by harmonic inversion of the periodic orbit signal C\{t) (see Eq. 
3.69), 

CiO = ~i E -l-J&W - U = i E djkf [kf - kf] e-*M , (3.71) 

po *-po j 



with dj the multiplicities of states. The signal C\(l) in Eq. 3.71 can be inverted 
as a single function as has been done in Ref. [86], where the accuracy of the 
eigenenergies was improved by one to several orders of magnitude, apart from 
the nearly degenerate states marked by arrows in Fig. 35. Here we go one step 



further and use Eq. 3.71 as part of a 3 x 3 cross-correlation signal. For the two 
other diagonal components of the cross-correlation matrix we use the identity 
operator, A\ = I, and the distance from the origin, A 2 = r. By applying the 
cross-correlation technique of Section 3.5 we obtain both the zeroth and first 
order h expansion of the eigenenergies and the diagonal matrix elements of the 
chosen operators simultaneously from one and the same harmonic inversion 
procedure. Furthermore, we can now even resolve the nearly degenerate states. 
The spectrum of the integrated differences of the density of states / Ag{k)dk 
obtained by harmonic inversion of the 3x3 cross- correlation matrix with sig- 
nal length £ max = 150 is shown in Fig. 36. The squares mark the spectrum 
for Ag(k) = g^(k) — g^(k) obtained from the harmonic inversion of the sig- 
nal Ci(s). For comparison the crosses present the same spectrum but for the 
difference Ag(k) = g cx (k) — g EBK (k) between the exact quantum mechanical 
and the EBK-spectrum. The deviations between the peak heights exhibit the 
contributions of terms of the h expansion series beyond the first order approx- 
imation. The peak heights of the levels in Fig. 36 (solid lines and crosses) are, 
up to a multiplicity factor for the degenerate states, the shifts Ak between 
the zeroth and first order semiclassical approximations to the eigenvalues k. 
The zeroth and first order eigenvalues, k^ and k^ = k^ + Ak are presented 
in Table 12 for the 40 lowest eigenstates. The zeroth order eigenvalues, k^°\ 
agree within the numerical accuracy with the results of the torus quantiza- 
tion, Eq. 2.71 (see eigenvalues k EBK in Table 12). However, the semiclassical 
eigenvalues deviate significantly, especially for states with low radial quantum 
numbers n, from the exact quantum mechanical eigenvalues k cx in Table 12. 
By contrast, the semiclassical error of the first order eigenvalues, A;'- 1 -', is by 
orders of magnitude reduced compared to the lowest order approximation. 

An appropriate measure for the accuracy of semiclassical eigenvalues is the 
deviation from the exact quantum eigenvalues in units of the average level 
spacings, (Ak) av = l/g(k). Fig. 37 presents the semiclassical error in units of 
the average level spacings (Ak) ay ~ 4/k for the zeroth order (diamonds) and 
first order (crosses) approximations to the eigenvalues. In the zeroth order 
approximation the semiclassical error for the low lying states is about 3 to 
10 percent of the mean level spacing. This error is reduced in the first order 
approximation by at least one order of magnitude for the least semiclassical 
states with radial quantum number n — 0. The accuracy of states with n > 1 
is improved by two or more orders of magnitude. 

Finally, we want to note that the small splittings between the nearly degen- 
erate states are extremely sensitive to the higher order h corrections. E.g., in 
the zeroth order approximation the splitting between the two states around 
k w 11.05 is A&(°) = 6 x 10~ 4 . In the first order approximation the splitting 
between the same states is AA^ 1 ) = 0.0242, which is very close to the exact 
splitting AA; CX = 0.0217. The accuracy obtained here for the first order approx- 
imations to the nearly degenerate states goes beyond the results presented in 



Ref . [86] , and is achieved by the combined application of the methods intro- 
duced in Sections 3.5 and 3.6, i.e., the harmonic inversion of cross-correlation 
functions and the analysis of the h expanded periodic orbit sum, respectively. 



3. 8 Semiclassical calculation of transition matrix elements for atoms in ex- 
ternal fields 



The interpretation of photoabsorption spectra of atoms in external fields is a 
fundamental problem of atomic physics. Although the "exact" quantum me- 
chanics accurately describes the energies and transition strengths of individual 
levels it has completely failed to present a simple physical picture of the long 
ranged modulations which have been observed in early low resolution spectra 
of barium atoms in a magnetic field [141] and later in the Fourier transform 
recurrence spectra of the magnetized hydrogen atom [142,71]. However, the 
long ranged modulations of the quantum photoabsorption spectra can be nat- 
urally interpreted in terms of the periods of classical closed orbits starting at 
and returning back to the nucleus where the initial state is localized. The link 
between the quantum spectra and classical trajectories is given by closed orbit 
theory [90-93] which describes the photoabsorption cross section as the sum of 
a smooth part and the superposition of sinusoidal modulations. The frequen- 
cies, amplitudes, and phases of the modulations are directly obtained from the 
quantities of the closed orbits. When the photoabsorption spectra are Fourier 
transformed or analyzed with a high resolution method (see Sections 2.1 to 
2.3) the sinusoidal modulations result in sharp peaks in the Fourier transform 
recurrence spectra, and closed orbit theory has been most successful to explain 
quantum mechanical recurrence spectra qualitatively and even quantitatively 
in terms of the closed orbits of the underlying classical system [46-48,51]. 

However, up to now practical applications of closed orbit theory have always 
been restricted to the semiclassical calculation of low resolution photoabsorp- 
tion spectra for the following two reasons. First, the closed orbit sum requires, 
in principle, the knowledge of all orbits up to infinite length, which are usu- 
ally not available from a numerical closed orbit search, and second, the infinite 
closed orbit sum suffers from fundamental convergence problems [90-93]. It is 
therefore commonly accepted that the calculation of individual transition ma- 
trix elements (<f)i\D\ijjf) of the dipole operator D, which describe the transition 
strengths from the initial state \<pi) to final states 1^/)) is a problem beyond 
the applicability of the semiclassical closed orbit theory, i.e., is the domain of 
quantum mechanical methods. 

In this section we disprove this common believe and demonstrate that indi- 
vidual eigenenergies and transition matrix elements can be directly extracted 
from the quantities of the classical closed orbits. To that end, we slightly gen- 
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eralize closed orbit theory to the semiclassical calculation of cross-correlated 
recurrence functions. We then adopt the cross-correlation approach introduced 
in Section 3.5 to harmonically invert the cross-correlated recurrence signal and 
to extract the semiclassical eigenenergies and transition matrix elements. Re- 
sults will be presented for the photo excitation of the hydrogen atom in a 
magnetic field. 

The oscillator strength / for the photo excitation of atoms in external fields 
can be written as 

/(£) = _ ^ ( £ _ Ei) Im (fclDG+Dlfc) ? (3.72) 

7T 

where \<pi) is the initial state at energy E^, D is the dipole operator, and 
G\ the retarded Green's function of the atomic system. The basic steps for 
the derivation of closed orbit theory are to replace the quantum mechanical 
Green's function in (3.72) with its semiclassical Van Vleck-Gutzwiller approx- 
imation and to carry out the overlap integrals with the initial state \<f>i). Here 
we go one step further by introducing a cross-correlation matrix 

g aa > = (<t>a\DG+D\<t> a ,) (3.73) 

with \<f) a ), ol — 1, 2, . . . , L a set of independent initial states. As will be shown 
below the use of cross-correlation matrices can considerably improve the con- 
vergence properties of the semiclassical procedure. In the following we will con- 
centrate on the hydrogen atom in a magnetic field with 7 = 5/(2.35 x 10 5 T) 
the magnetic field strength in atomic units. As discussed in Section 2.1 the 
system has a scaling property, i.e., the shape of periodic orbits does not de- 
pend on the scaling parameter, w = 7 -1 / 3 = h^, and the classical action 
scales as S = sw with s the scaled action. As, e.g., in Ref. [47] we consider 
scaled photoabsorption spectra at constant scaled energy E = E^~ 2 I Z as a 
function of the scaling parameter w. We choose dipole transitions between 
states with magnetic quantum number m — 0. Note that the following ideas 
can be applied in an analogous way to atoms in electric fields. Following the 
derivation of Refs. [91,93] the semiclassical approximation to the fluctuating 
part of g aa i in Eq. 3.73 reads 



w 



■(2 



7T 



,5/2 



CCt 

mi2\ 



-1/2 j^ v ; Jsin^sin^ 



7 

x y a (^T)ya'(^f) e l ^° w -^° + V , (3.74) 

with s co and /i co the scaled action and Maslov index of the closed orbit 
(co), m c ° 2 an element of the monodromy matrix, and df and i9SP the ini- 
tial and final angle of the trajectory with respect to the magnetic field axis. 
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The angular functions y a ($) depend on the states \<f> a ) and the dipole op- 
erator D and are given as a linear superposition of Legendre polynomials, 
X*C#) — J2iBiaPi(cos$) with usually only few nonzero coefficients E\ a with 
low I. Explicit formulae for the calculation of the coefficients can be found in 
Refs. [91,93] and in Appendix B. The problem is now to extract the semiclas- 
sical eigenenergies and transition matrix elements from Eq. 3.74 because the 
closed orbit sum does not converge. The Fourier transformation of w 1 ^ 2 g^ a ,(w) 
yields the cross-correlated recurrence signals 

C b :As) = T,^Z^(s-s co ), (3.75) 

CO 

with the amplitudes 

A c a ° a > = ~ \— — A/sin ■&? sin 0f ^(0X^7) e<~^ 0+ ^ (3.76) 

'i „co| V J J 

'12 1 



being determined exclusively by closed orbit quantities. The corresponding 
quantum mechanical cross-correlated recurrence functions, i.e., the Fourier 
transforms of w 1 ^ 2 g'^i{w) read 

CZ(s) = -i52b ak ba' k e- iWkS , (3.77) 

k 



with Wk the eigenvalues of the scaling parameter, and 

b ak = wl'\<t> a \D\i) k ) (3.78) 

proportional to the transition matrix element for the transition from the initial 
state \4> a ) to the final state \ipk)- 

The method to adjust Eq. (3.75) to the functional form of Eq. (3.77) for a set of 
initial states \<f> a ), a — 1, 2, . . . , L is the harmonic inversion of cross-correlation 
functions as discussed in Section 3.5 and Appendix A. 2. We now demonstrate 
the method of harmonic inversion of the cross-correlated closed orbit recur- 
rence functions (3.75) for the example of the hydrogen atom in a magnetic 
field at constant scaled energy E = —0.7. This energy was also chosen for 
detailed experimental investigations on the helium atom [48]. We investigate 
dipole transitions from the initial state |0i) = |2p0) with light polarized par- 
allel to the magnetic field axis to final states with magnetic quantum number 
m — 0. For this transition the angular function in Eq. 3.76 reads (see Ap- 
pendix B) y^-d) = (27r)- 1 / 2 2 7 e- 4 (4cos 2 ^- 1). For the construction of a 2 x 2 
cross-correlated recurrence signal we use for simplicity as a second transition 
formally an outgoing s-wave, i.e., D\(p 2 ) oc Y 0t0 , and, thus, 3M"#) — const. 
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A numerical closed orbit search yields 1395 primitive closed orbits (2397 or- 
bits including repetitions) with scaled action s/2n < 100. With the closed 
orbit quantities at hand it is straightforward to calculate the cross-correlated 
recurrence functions in (3.75). The real and imaginary parts of the complex 
functions Cf[(s), C^(s), and C^(s) with s/2it < 50 are presented in Figs. 38 
and 39, respectively. Note that for symmetry reasons C|i(s) = Cf^is). 

We have inverted the 2x2 cross-correlated recurrence functions in the region 
< s/2tt < 100. The resulting semiclassical photoabsorption spectrum is 
compared with the exact quantum spectrum in Fig. 40a for the region 16 < 
w < 21 and in Fig. 40b for the region 34 < w < 40. The upper and lower parts 
in Fig. 40 show the exact quantum spectrum and the semiclassical spectrum, 
respectively. Note that the region of the spectrum presented in Fig. 40b belongs 
well to the experimentally accessible regime with laboratory field strengths 
B = 6.0 T to B = 3.7 T. The overall agreement between the quantum and 
semiclassical spectrum is impressive, even though a line by line comparison 
still reveals small differences for a few matrix elements. It is important to note 
that the high quality of the semiclassical spectrum could only be achieved 
by our application of the cross-correlation approach. For example, the two 
nearly degenerate states at w — 36.969 and w = 36.982 cannot be resolved 
and the very weak transition at w — 38.894 with (2pO\D\ipf) 2 = 0.028 is not 
detected with a single (lxl) recurrence signal of the same length. However, 
these hardly visible details are indeed present in the semiclassical spectrum 
in Fig. 40b obtained from the harmonic inversion of the 2x2 cross-correlated 
recurrence functions. 

To summarize, we have demonstrated that closed orbit theory is not restricted 
to describe long ranged modulations in quantum mechanical photoabsorption 
spectra of atoms in external fields but can well be applied to extract individual 
eigenenergies and transition matrix elements from the closed orbit quantities. 
This is achieved by the high resolution spectral analysis (harmonic inversion) 
of cross-correlated closed orbit recurrence signals. For the hydrogen atom in a 
magnetic field we have obtained individual transition matrix elements between 
low lying and highly excited Rydberg states solely from the classical closed 
orbit data. 



4 Conclusion 



In summary, we have shown that harmonic inversion is a powerful tool for the 
analysis of quantum spectra, and is the foundation for a novel and universal 
method for periodic orbit quantization of both regular and chaotic systems. 

The high resolution analysis of finite range quantum spectra allows to circum- 



vent the restrictions imposed by the uncertainty principle of the conventional 
Fourier transformation. Therefore physical phenomena can directly be revealed 
in the quantum spectra which previously were unattainable. Topical examples 
are the study of quantum manifestations of periodic orbit bifurcations and 
catastrophe theory, and the uncovering of symmetry breaking effects. The 
investigation of these phenomena provides a deeper understanding of the rela- 
tion between quantum mechanics and the dynamics of the underlying classical 
system. The high resolution technique is demonstrated in this work for nu- 
merically calculated quantum spectra of, e.g., the hydrogen atom in external 
fields, three disk pinball scattering, and the circle billiard. Theoretical spectra 
are especially suited for the high resolution analysis because of the very high 
accuracy of most quantum computational methods. In principle, harmonic in- 
version may be applied to experimental spectra as well, e.g., to study atoms in 
external fields measured with the technique of scaled- energy spectroscopy [45- 
55]. However, the exact requirements on the precision of experimental data 
to achieve high resolution recurrence spectra beyond the limitations of the 
uncertainty principle are not yet known. It will certainly be a challenge for fu- 
ture experimental work to verify the quantum manifestations of bifurcations, 
which have been extracted here from theoretically computed spectra, using 
real systems in the laboratory. 

We have also introduced harmonic inversion as a new and general tool for 
semiclassical periodic orbit quantization. Here we briefly recall the highlights 
of our technique. The method requires the complete set of periodic orbits up 
to a given maximum period as input but does not depend on special prop- 
erties of the orbits, as, e.g., the existence of a symbolic code or a functional 
equation. The universality and wide applicability has been demonstrated by 
applying it to systems with completely different properties, namely the ze- 
ros of the Riemann zeta function, the three disk scattering problem, and the 
circle billiard. These systems have been treated before by separate efficient 
methods, which, however, are restricted to bound ergodic systems, systems 
with a complete symbolic dynamics, or integrable systems. The harmonic in- 
version technique allows to solve all these problems with one and the same 
method. The method has furthermore been successfully applied to the hydro- 
gen atom in a magnetic field as a prototype example of a system with mixed 
regular-chaotic dynamics. The efficiency of the method can be improved if 
additional semiclassical information obtained from a set of linearly indepen- 
dent observables is used to construct a cross-correlated periodic orbit sum, 
which can then be inverted with a generalized harmonic inversion technique. 
The cross-correlated periodic orbit sum allows the calculation of semiclassical 
eigenenergies from a significantly reduced set of orbits. Eigenenergies beyond 
the Gutzwiller and Berry- Tabor approximation are obtained by the harmonic 
inversion of the h expansion of the periodic orbit sum. When applied, e.g., 
to the circle billiard the semiclassical accuracy is improved by at least one to 
several orders of magnitude. The combination of closed orbit theory with the 



cross-correlation approach and the harmonic inversion technique also allows 
the semiclassical calculation of individual quantum transition strengths for 
atoms in external fields. 

Periodic orbit quantization by harmonic inversion has been applied in this 
work to systems with scaling properties, i.e., systems where the classical ac- 
tions of periodic orbits depend linearly on a scaling parameter, w. However, 
the method can even be used for the semiclassical quantization of systems 
with non-homogeneous potentials such as the potential surfaces of molecules. 
The basic idea is to introduce a generalized scaling technique with the inverse 
Planck constant w = 1/h as the new formal scaling parameter. The general- 
ized scaling technique can be applied, e.g., to the analysis of the rovibrational 
dynamics of the HO2 molecule [94]. By varying the energy of the system, 
the harmonic inversion method yields the semiclassical eigenenergies in the 
(E, w) plane. For non-scaling systems the semiclassical spectra can then be 
compared along the line with the true physical Planck constant, w = 1/h = 1, 
with experimental measurements in the laboratory. 
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A Harmonic inversion by filter-diagonalization 



In the following we give details about the numerical method of harmonic 
inversion by filter-diagonalization. We begin with the harmonic inversion of 
a single function and then extend the method to the harmonic inversion of 
cross-correlation functions. 



A.l Harmonic inversion of a single function 



The harmonic inversion problem can be formulated as a nonlinear fit (see, e.g., 
Ref. [67]) of the signal C(s) defined on an equidistant grid, 



C(nr) = J2 d k e~ inrWk , n = 0, 1, 2, ... TV, (A.l) 



with the set of generally complex variational parameters {w k ,dk}- (In this 
context the Discrete Fourier Transform scheme would correspond to a linear fit 
with iV amplitudes d k and fixed real frequencies w k = 2iik/NT, k = 1, 2, . . . N. 
The latter implies the so called "uncertainty principle", i.e., the resolution, 
defined by the Fourier grid spacing, Aw, is inversely proportional to the length, 
Smax = Nt, of the signal C(s).) The "high resolution" property associated with 
Eq. A.l is due to the fact that there is no restriction for the closeness of the 
frequencies w k as they are variational parameters. In Ref. [68] it was shown 
how this nonlinear fitting problem can be recast as a linear algebraic one 
using the filter-diagonalization procedure. The essential idea is to associate 
the signal c n with an autocorrelation function of a suitable dynamical system, 

%,U n %), (A.2) 



where ( • , • ) defines a complex symmetric inner product (i.e., no complex 
conjugation). The evolution operator can be defined implicitly by 

U = e~ irh = J2e- lTWk \r k )(T k \ , (A.3) 

fc=i 



where the set of eigenvectors {T^} is associated with an arbitrary orthonor- 
mal basis set and the eigenvalues of U are u k = e~ lTWk (or equivalently the 
eigenvalues of Cl are w k ). Inserting Eq. A.3 into Eq. A.2 we obtain Eq. A.l 
with 

4 = (T fc ,$ ) 2 , (A.4) 



which also implicitly defines the "initial state" $o- This construction estab- 
lishes an equivalence between the problem of extracting spectral information 
from the signal with the one of diagonalizing the evolution operator U = e~ tTn 
(or the Hamiltonian Cl) of the fictitious underlying dynamical system. The 
filter-diagonalization method is then used for extracting the eigenvalues of the 
Hamiltonian Cl in any chosen small energy window. Operationally this is done 
by solving a small generalized eigenvalue problem whose eigenvalues yield the 
frequencies in a chosen window. The knowledge of the operator Cl itself is 
not required, as for a properly chosen basis the matrix elements of Cl can 
be expressed only in terms of c n . The advantage of the filter-diagonalization 
procedure is its numerical stability with respect to both the length and com- 
plexity (the number and density of the contributing frequencies) of the sig- 
nal. Here we apply the method of Refs. [69,70] which is an improvement 
of the filter-diagonalization method of Ref. [68] in that it allows to signifi- 
cantly reduce the required length of the signal by implementing a different 
Fourier-type basis with an efficient rectangular filter. Such a basis is defined 
by choosing a small set of values (pj in the frequency interval of interest, 
Tw m in < ipj < rw max , j = 1,2,..., J, and the maximum order, M, of the 
Krylov vectors, $„ = e -m <&0j used in the Fourier series, 



M M 

n=0 n=0 



ttj = Vfa) = J2 e in<Pi ®n = E e m( ^- r % . (A.5) 



It is convenient to introduce the notations, 

U$ = U^( Vv w) = (tf fo), e-^%(yv)) , (A.6) 

for the matrix elements of the operator e~ iprn , and U' p % for the corresponding 
small J x J complex symmetric matrix. As such U' 1 - 1 denotes the matrix rep- 
resentation of the operator U itself and U' -* , the overlap matrix with elements 
(ty(<Pj),ty(<Pj>)), which is required as the vectors ^(<fj) are not generally or- 
thonormal. Now using these definitions we can set up a generalized eigenvalue 
problem, 

U (p) B fe = e- ipTWk XJ (0) B k , (A.7) 

for the eigenvalues e~ ipTWk of the operator e~ ip . The column vectors B& with 
elements Bj^ define the eigenvectors Y& in terms of the basis functions \Pj as 

T^E B A ( A -8) 

i=i 

assuming that the Size's form a locally complete basis. 
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The matrix elements (A. 6) can be expressed in terms of the signal c n , the 
explicit knowledge of the auxiliary objects Q, Y& or $ is n °t needed. Indeed, 
insertion of Eq. A. 5 into Eq. A. 6, use of the symmetry property, (\l/, £/$) = 
(Uty, $), and the definition of c n , Eq. A. 2, gives after some arithmetics 



M 



U (p \<p, <p') = (e-* - e'Vy 1 e -^ ^ e inip ' c n 



(A.9) 



n=0 



M 



-V^ 



e xwp c 



n+p 



JMtfi 



n=0 



2M 

E • 

n=M+l 



3 i(ra-M-lW 



^n+p 



+ e 

2M 

E 

n=0 



iMip' 



2M 

E 



i{n-M-\)if> 



Cn+z 



cp^cp' 



n=M+l 

U^(ip, ip) = £(M - \M - n\ + l)e m ^c n+p . 



(Note that the evaluation of U^ requires the knowledge of c n for n = p, p + 
l,...,N = 2M + p.) 

The generalized eigenvalue problem (A. 7) can be solved by a singular value 
decomposition of the matrix U^ -* , or more accurately by application of the QZ 
algorithm [143], which is implemented, e.g., in the NAG library. Each value of 
p yields a set of frequencies w k and, due to Eqs. A. 4, A. 5 and A. 8, amplitudes, 



J M 

J2B jk J2c n e^ 

Kj = l n=0 



(A.10) 



Note that Eq. A. 10 is a functional of the half signal c n ,n = 1,2,..., M. An 
even better expression for the coefficients d k (see Ref. [70]) reads 



di 



M + 



1 J 



1 J 



M + 



-Y,B jk UW{tpj,Wk) 

1 3=1 



(A.11) 



with U^((fj,Wk) defined by Eq. A.9. Eq. A. 11 is a functional of the whole 
available signal c n , n — 0, 1, . . . , 1M and therefore sometimes provides more 
precise results than Eq. A. 10. 

The converged w k and d k should not depend on p. This condition allows us to 
identify spurious or non-converged frequencies by comparing the results with 
different values of p (e.g., with p — 1 and p = 2). We can define the simplest 
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error estimate e as the difference between the frequencies Wk obtained from 
diagonalizations with p — 1 and p = 2, i.e. 

e = \w% =1) -w$ =2) \. (A.12) 

A. 2 Harmonic inversion of cross-correlation functions 



We now consider a cross-correlation signal, i.e., a D x D matrix of signals 
defined on an equidistant grid (a, a' = 1,2, ... ,D): 



c nn >(n 



C aa ,(nT) = J2ba,kb a i, k e- inTWk , n = 0, 1, 2, . . . , N . (A.13) 



[We choose r = 1 for simplicity in what follows.] Each component of the 
signal c aa i{n) contains the same set of frequencies Wk, and the amplitudes 
belonging to each frequency are correlated, i.e., d aoi i^ = b a ^b a ',k with only D 
(instead of D 2 ) independent parameters b a ^. As for the harmonic inversion of 
a single function the nonlinear problem of adjusting the parameters {wk, & a ,fc} 
can be recast as a linear algebra one using the filter-diagonalization procedure 
[68,137,138]. The cross-correlation signal (A.13) is associated with the cross- 
correlation function of a suitable dynamical system, 



($ a ,£T<lv) , (A.14) 



with the same complex symmetric inner product as in Eq. A. 2, and the evo- 
lution operator U defined implicitly by Eq. A. 3. The extension from Eq. A. 2 
to Eq. A.14 is that the autocorrelation function (§ ,U n § ) built of a single 
state $o is replaced with the cross-correlation function ($ a , U n Q a r) built of 
a set of D different states $ Q . Inserting (A. 3) into Eq. (A.14) we obtain Eq. 
(A.13) with 



Ja,k 



(T fc ,$ a ), (A.15) 



which now implicitly defines the states $ Q . After choosing a basis set in anal- 
ogy to Eq. A. 5, 



M 



* QJ = \I> Q (^) = £ e in to- n >$ Q , (A.16) 



n=0 



and introducing the notations 

U (p) ee U^ a , f = U®(<p s , tpy) = (* a ( W ), e-^Mv(^O) , (A.17) 
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for the small matrix of the operator e tpQ in the basis set (A. 16), we can set 
up a generalized eigenvalue problem, 

U (p) B fe = e - ipWk U {0) B k , 



for the eigenvalues e~ tpWk of the operator e~ ipn , which is formally identical with 
Eq. A. 7. The matrix elements (A. 17) can be expressed in terms of the signal 
Caa',n- The following expression for the matrix elements of U^ is derived in 
complete analogy with Refs. [69,70] with the additional indices a, a' being the 
only difference, 



M 
r(p) 



U™(<p, <//) = (e-* - e-^'y 1 [e-* £ e in *'c aa ,(n + p) 

M 

-e-^'^e^Caain + p) 



n=0 

2M 

_ e iM V j- e^-M-^'c^in + p) 

n=M+l 
2M 

+e i Mv > J2 e^-^cwln + p)] ,<P^<p', (A- 18) 

n=M+l 
2M 

[/£, (y?, </?) = £(M - |M - n| + l)e^w(n + p) . (A.19) 

ra=0 

Given the cross-correlation signal c aa >(n), the solution of the generalized eigen- 
value problem (A. 7) yields the eigenfrequencies Wk and the eigenvectors B&. 
The latter can be used to compute the amplitudes, 



'ak 



' "iwEE B^U^i^Wk + n) , (A.20) 



1 _ e -(M+l) 7 t 

j=l a'=l 



where the adjusting parameter 7 is chosen so that U^((fj, Wk + ij) is numer- 
ically stable [138]. One correct choice is 7 = — Im Wk for Im Wk < and 7 = 
for Im Wk > 0. 



B Angular function y m {$) 



In Section 2.4.2 we have presented Eq. 2.36 as the final result of closed orbit 
theory for the semiclassical photoabsorption spectrum of the hydrogen atom 
in a magnetic field. Here we define explicitly the angular function 3-m($) m 
Eq. 2.36. 
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The angular function !V m (i?) solely depends on the initial state ipi and the 
dipole operator D and is a linear superposition of spherical harmonics: 

oo 

y m (#) = E (-i)% m Y tm (&,o) ■ (b.i) 

i=\m\ 



The coefficients Be m are defined by the overlap integrals 

B tm = Jd 3 x(D^)(x)^2/r~J 2e+1 (v^) Y e * m (#, tp) (B.2) 



(with J„(a;) the Bessel functions) and can be calculated analytically [93]. For 
excitations of the ground state ipi = |lsO) with 7r-polarized light (i.e. dipole 
operator D = z) the explicit result is 

y {ff) = -ir- 1/2 2 3 e- 2 cos ■& , (B.3) 



and for ipi = \2p0), i.e., the initial state in many spectroscopic measurements 
on hydrogen [45-47] we obtain 

y (&) = (2ny 1 / 2 2 7 e~ 4 {a cos 2 d - l) . (B.4) 



C Catastrophe diffraction integrals 



Here we give some technical details about the numerical calculation and the 
asymptotic expansion of catastrophe diffraction integrals for the hyperbolic 
umbilic and the butterfly catastrophe. 



C.l The hyperbolic umbilic catastrophe 

The catastrophe diffraction integral of the hyperbolic umbilic (Eq. 2.28) reads 
V(x,y) = I dp ( dqe l[p3W+y{p+q)2+x{p+q)] . (C.l) 



-co — oc 



By substituting 



y 
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q = s 2 



we obtain 



y(x,y) = e l Ml y2 - x )h(x-y 2 ,2y) , 



(C.2) 



with 



-oo +00 



$(£,77) = I ds x f ds 2 e i ^ +s ^ + ^ sl+S2)+r ' slS2 



-00 — 00 



(C.3) 



The integral $(£, 77) can be expanded into a Taylor series around rj = 0. Using 



CT$(£,77) 



-00 +00 



di] r ' 



??=0 



dsi / ds 2 i n (sis 2 )V [s ? +s 2+«( s i+ s 2)] 



-co —CO 

+ CO 



i n ds s V (s3+?s) 



(C.4) 



and solving the one-dimensional integrals [144] 



+00 



dss n e i{ S t+ts) 



l^(-fl* / fc + n + T 
3fro *' V 3 t 



i(n+l-2k)ir/6 , /-i \n -i(n+l— 2k)ir/6 //~i c'i 



we finally obtain 



•fco-ss^g^p^: 



i(n+l-2k)-tt/S , /_ i \n -i(n+l-2fc)ir/6' 



(C.6) 



which is a convergent series for all £ and 77. 



C.2 TTie butterfly catastrophe 



The uniform phase integral \l/(x, y) of the butterfly catastrophe is expanded 
in a two-parametric Taylor series around x = y = 0: 
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*(x, y) = / exp -i(xf + yt 4 + *' 



,4 i ^ 



r/t 



-oo 

oo oo 



+00 



EE^^T /t— exp 



n=0 m=0 



i n+m n\ m\ 



-if 



dt . 



With the substitution z = t 2n+4m+1 we obtain [144] 



(C.7) 



+00 



£ 2n+4m exp 



-if 



dt 



2n + 4m + 1 



exp 



• 6/(2n+4m+l) 



= -exp 
and finally 



2n + 4m + 1 ' 

-l 7T 

12 



fl.2 

2n + 4m + 1 
6 



^(^,l/) = 3 ex P 



7T 



y, y, 1 r ^2n + 4m + l 



12J^ ^o n!m! ^ 

2 l\ n / r 5 

y exp 



6 



-Z— 7T 
3 



x I x exp 
which is a convergent series for all x and ?/. 



2-7T 

6 . 



(C.8) 



The asymptotic behavior of ^(x, y) for x — ► ±00 is obtained in a stationary 
phase approximation to Eq. (C.7) with the stationary points t being defined 
by 



t (6f + 4yt 2 + 2x) = 0. 



(C.9) 



a) x — ► —00: 

There are three real stationary points given by 

2 1 /ft r~ 

t = and t = --y + \Jgy 3^ 



and we obtain 



V(x,y) 1 



exp 



— ixt z 



dt + 2 exp -i(xtl + yf\ + f\ 
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x / exp 

— oo 



4\ .2 



-i(x + 6yt 2 + 15Q t 



dt 



7[ 



—X 



x exp < i 



exp 



,7T 



+ 



7T 



\ -* + iv 2 - y^JIv 2 - i 



X 



■ , 1 2 1 \ 3/2 1 /2 2 

2-1/ x v [-V — x 



(CIO) 



b) x — *■ +oo: 

The only real solution of (C.9) is t = and \&(:r,i/) is approximated as 



+oo 



*(x,y) 



exp 



-ixt z 



fit 



TV 

— exp 



7T 

-I — 

4 



(Cll) 



In the case of y < and a; ~ |i/ 2 the complex zeros 



tn = =b 



-nV±i\h; x 



1 1 



\"3 



of Eq. (C.9) are situated close to the real axis. When one considers these 
complex zeros in a stationary phase approximation Eq. (Cll) is modified by 
an additional term exponentially damped for large x: 
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Fig. 1. (a) Experimental recurrence spectra of the hydrogen atom in a magnetic 
field in overlay form. Even-parity, m = final state, (b) Scaled action of closed 
classical orbits through origin. (From Ref. [45].) 
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Fig. 2. Closed classical orbits of the hydrogen atom in a magnetic field in (p, z) 
projection, (a) Vibrator, (b) rotator, and (c) "exotic" orbits. (From Ref. [46].) 
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Fig. 3. Experimental recurrence spectra of lithium in an electric field in overlay 
form. The curves mark the scaled action of the parallel orbit and its repetitions as 
a function of scaled energy. The open circles mark the bifurcation points of closed 
orbits, whose shapes are shown along the bottom. (From Ref. [52].) 
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Fig. 4. Recurrence spectrum (Fourier transform) for the density of states of the 
hydrogen atom in a magnetic field at scaled energy E = E"/~ 2 ' 3 = —0.1. 
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Fig. 5. (a) Recurrence spectrum for the density of states of the hydrogen atom in 
a magnetic field at scaled energy E = E^~ 2 '^ = —0.1. Smooth line: Conventional 
Fourier transform. Solid stick spectrum: High resolution quantum recurrence spec- 
trum. Dashed sticks (hardly visible under solid sticks): Recurrence spectrum from 
semiclassical periodic orbit theory, (b) Complex actions. Crosses and squares are the 
quantum and classical results, respectively. 
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Fig. 6. Recurrence spectrum for the density of states of the hydrogen atom in a 
magnetic field at scaled energy E = Ej~ 2 ' 3 = —0.1. Smooth line: Conventional 
Fourier transform. Solid stick spectrum and crosses: High resolution quantum re- 
currence spectrum. Dashed sticks (hardly visible under solid sticks) and squares: 
Recurrence spectrum from semiclassical periodic orbit theory. The recurrence peaks 
are identified by periodic orbits drawn in semiparabolical coordinates n = (r+z) 1 ' 2 , 
v={r- z) 1 / 2 . (From Ref. [79].) 
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Fig. 7. Schematic view of the four-disk scattering problem. (From Ref. [102]. 
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Fig. 8. Difference As between the classical action of the four periodic orbits involved 
in the bifurcations. Dashed lines: Analytical fits related to the hyperbolic umbilic 

catastrophe. Inset: Graphs of periodic orbits 0-+ — and its time reversal + (solid 

line), 00+- (dashed line), and +++ (dash-dotted line) drawn in semiparabolical 

coordinates \x = (r + z) l < 2 , v = (r — z) 1 ' 2 . (From Ref. [44].) 
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Same as Fig. 8 but for the determinant det(M — /) of the periodic orbits. 
Ref. [44].) 
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Fig. 10. Semiclassical amplitudes (absolute values) for magnetic field strength (a) 
7 = 10 , (b) 7 = 10 -8 , and (c) 7 = 10 -9 in units of the time period To. Dashed line: 
Amplitudes of the standard semiclassical trace formula. Dash-dotted line: Ghost or- 
bit contribution. Solid line: Uniform approximation of the hyperbolic umbilic catas- 
trophe. (From Ref. [44].) 
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Fig. 11. High resolution recurrence spectra at scaled energy E = —0.1. Solid peaks: 
Part of the quantum recurrence spectra. Dash-dotted peak: Classical ghost orbit 
contribution. Dashed peaks: Uniform approximation of the hyperbolic umbilic catas- 
trophe. (From Ref. [44].) 
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Fig. 12. (a) Real and imaginary part of starting angle #j for closed orbits related 
to the bifurcating scenario of the (period doubled) perpendicular orbit i?2- (b) Dif- 
ference As/2tt between the classical action of the (period doubled) perpendicular 
orbit i?2 an d real and ghost orbits bifurcating from it. (c) Monodromy matrix ele- 
ment mi2 of the perpendicular orbit R2 and orbits bifurcating from it. Dashed lines: 
Analytical fits (see text). (From Ref. [37].) 
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Fig. 13. (a) Butterfly catastrophe of returning orbits (in rotated semiparabolic coor- 
dinates) related to the bifurcation of the perpendicular orbit at E\, = —0.31618537. 
(b) Magnification of the marked region close to the nucleus. There are one, three, 
or five orbits returning to each point (£,77) in coordinate space. (From Ref. [37].) 
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Fig. 14. (a) Recurrence spectrum for the photoabsorption cross section of the hy- 
drogen atom in a magnetic field at scaled energy E = E^ -2 '^ = —0.35. Transition 
|2p0) — ► \m 7rz = + ). Smooth line: Conventional Fourier transform. Solid stick 
spectrum and crosses: High resolution quantum recurrence spectrum. Dashed sticks 
and squares: Recurrence spectrum from semiclassical closed orbit theory. The two 
strongest recurrence peaks are identified by a real and complex ghost orbit which 
are presented as insets. The solid and dashed lines in the insets are the real and 
imaginary part in semiparabolical coordinates fi = (r + z) 1 ' 2 , v = (r — z) 1 ' 2 . (b) 
Complex actions. Crosses and squares are the quantum and classical results, respec- 
tively. (FromRef. [79].) 
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Fig. 15. Segments of (a) experimental and (b) theoretical scaled action spectra for 
the hydrogen atom in crossed magnetic and electric fields at constant scaled energy 
E = -0.15. (From Ref. [54].) 
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Fig. 16. (a) Classical trajectories for the hydrogen atom in crossed fields in the 
plane perpendicular to the magnetic field axis. Scaled energy E = —0.15 and field 
strength / = 0.012. (b) Closeup of the returning part of the trajectories as they 
approach the nucleus. (From Ref. [54].) 
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Fig. 17. Recurrence spectrum for the photo-absorption cross section of the hydrogen 
atom in crossed magnetic and electric fields at scaled energy E = E"/~ 2 ' 3 = —0.5 
and field strength / = F^ -4 / 3 = 0.02. Transition |2p0) — ► \7r z = +1). Smooth line: 
Conventional Fourier transform. Solid stick spectrum and crosses: High resolution 
quantum recurrence spectrum. Dashed sticks and squares: Recurrence spectrum 
from semiclassical closed orbit theory. The two recurrence peaks around s/2ir = 0.87 
are identified by closed orbits in the (x, y) plane presented as insets in the figure. 
(FromRef. [79].) 
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Fig. 18. Part of the difference spectrum Ag(k) = g Cim (k) — g sc (k) between the 
quantum and the semiclassical density of states for the circle billiard with radius 
R = l. 
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Fig. 19. First order ft correction terms of the semiclassical periodic orbit sum for the 
circle billiard. Solid line: Semiclassical Theory. Crosses: Harmonic inversion analysis 
of the difference spectrum k 1 ' 2 [Q qm (k) — g sc (k)]. The periodic orbits are marked by 
the numbers (m^,,?rt r ). 
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(b) 



Fig. 20. The scattering geometry for the three-disk system, (a) The three disks with 
12, 123, and 121313232 cycles indicated, (b) The fundamental domain, i.e., a wedge 
consisting of a section of a disk, two segments of symmetry axis acting as straight 
mirror walls, and an escape gap. The above cycles restricted to the fundamental 
domain are now the 0, 1, and 100 cycle. (From Ref. [9].) 
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Fig. 21. Spectra of the three disk scattering system {A\ subspace) with radius 
R = 1 and distance d = 6 between the disks, (a) Fluctuating part of the quantum 
mechanical density of states, (b) Difference spectrum k[g qm (k) — £ SC (A;)] between 
the quantum and semiclassical density of states, (c) Quantum (+) and semiclassical 
(x) A\ resonances of the Green's function. 
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Fig. 22. "Recurrence" function —iC a {s) for the Riemann zeros which has been 
analyzed by harmonic inversion, (a) Range < s < 6, (b) short range around 
s = 13. The (5-functions have been convoluted by a Gaussian function with width 
a = 0.0003. Dashed line: Smooth background C{s) = ie s < 2 resulting from the pole 
of the zeta function at w = i/2. (From Ref. 
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Fig. 23. Estimated number of converged zeros of the Riemann zeta function, which 
can be obtained by harmonic inversion for given number of primes p. (From Ref. 

[83].) 
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Fig. 24. The spectral staircase M + (E) and its semiclassical approximation J\f+(E) 
for the hyperbolic billiard. N£(E) was calculated using 101265 periodic orbits. 
(FromRef. [65].) 
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Fig. 25. The function cos{-7rA/' s ^(-E')} for the hyperbolic billiard. N£(E) was evalu- 
ated as in Fig. 24. The triangles mark the positions of the true quantum mechanical 
energies. (From Ref. [65].) 
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Fig. 26. The cycle-expanded functional determinant for the m n = + subspace of 
the anisotropic Kepler problem. The quantum eigenvalues are marked with vertical 
bars on the real axis, (a) All orbits up to length 4 (8 in number) included, (b) All 
orbits up to length 8 (71) included. (From Ref. [114].) 
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Fig. 27. The cycle-expanded functional determinant for the closed three disk billiard 
(disk radius R = 1). The vertical bars mark the exact quantum eigenvalues, (a) All 
orbits up to length 2 (3 in number) included, (b) All orbits up to length 3 (5) 
included. (From Ref. [114].) 
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Fig. 28. Three disk scattering system (A\ subspace) with R = 1, d = 6. (a) Periodic 
orbit recurrence function, C{L). The signal has been convoluted with a Gaussian 
function of width a = 0.0015. (b) Semiclassical resonances. (From Ref. [82].) 
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Fig. 29. Three disk scattering system (A\ subspace) with R = 1, d = 6. Semiclassical 
resonances obtained by harmonic inversion of the periodic orbit recurrence function 
C(L) with L < 24 (squares) and L < 50 (crosses). 
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Fig. 30. Three disk scattering system {A\ subspace) with R = 1, d = 2.5. (a) 
Periodic orbit recurrence function, C(L). The signal has been convoluted with a 
Gaussian function of width a = 0.0003. (b) Semiclassical resonances. (From Ref. 

[82].) 
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Fig. 31. Smoothed fluctuating part of the density of states for the hydrogen atom 
in a magnetic field at scaled energy E = —0.2. The quantum results (thick lines) 
are smoothed over the first 100 eigenstates. Semiclassical results (thin lines) are 
obtained by including (a) two orbits (three contributions including repetitions) 
with scaled actions s/2ir < 1.33, and (b) 13 orbits (19 including repetitions) with 
s/2tt < 3. The lowest quantum eigenvalues are marked as vertical bars. (From Ref. 
[78].) 
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Fig. 32. (a) Action s/2ttM2 and (b) second derivative g" E as a function of the 
frequency ratio MijM\ for the rotator (solid lines) and vibrator (dashed lines) type 
motion of the hydrogen atom in a magnetic field at scaled energy E = —0.4. Inset: 
Poincare surface of section in semiparabolical coordinates (/U,p M ; v = 0). (From Ref. 

[84].) 
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Fig. 33. Semiclassical recurrence signal C sc (s) for the hydrogen atom in a magnetic 
field at scaled energy E = —0.4. Solid and dashed sticks: Signal from the rotator 
and vibrator type motion, respectively. (From Ref. [84].) 
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Fig. 34. (a) Semiclassical and quantum mechanical spectrum of the hydrogen atom 
in a magnetic field at scaled energy E = —0.4. (b) Distribution P(d) of the semi- 
classical error in units of the mean level spacing, d = (u> qm — w S c)/Aw a v for the 
lowest 106 eigenvalues. (From Ref. [84].) 
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Fig. 36. Integrated difference of the density of states, / Ag(k)dk, for the cir- 
cle billiard with radius R = 1. Crosses: Ag(k) = g cx (k) — g (k). Squares: 
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Fig. 37. Semiclassical error |/c' ' — A; cx | (diamonds) and |fc") — k cx \ (crosses) in units 
of the average level spacing (A/c) av ss 4/fc. 



138 



co 



O 

CD 
DC 



200 

100 

-100 
-200 h 



~\ r 





200 


3* 
O 


100 



CD 
DC 


-100 



-200 





200 


3- 

o 


100 



CD 
DC 


-100 



-200 



"i 1 r 



j i i L 



~\ 1 r 



J L 



~\ r 



\\\ 



! 't 



iwl 



p w 



J I I I I L 



(a) 



ft. 



J I I L 



"i 1 1 1 1 1 1 r 



(b) 



"i 1 1 1 1 1 r 



(C) 



^ 



TPfl 



5 10 15 20 25 30 35 40 45 50 

S/27U 

Fig. 38. Real parts of the cross-correlated recurrence functions for the hydrogen 
atom in a magnetic field at constant scaled energy E = —0.7. (From Ref. [89].) 
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Fig. 39. Same as Fig. 38 but for the imaginary parts of the cross-correlated recur- 
rence functions. 
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Fig. 40. Quantum (upper part) and semiclassical (lower part) photoabsorption spec- 
tra of the hydrogen atom in a magnetic field at scaled energy E = —0.7. Transition 
matrix elements (2pO\D\ipf) 2 for dipole transitions with light polarized parallel to 
the magnetic field axis. (From Ref. [89].) 
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6 Tables 



Table 1 

Hydrogen atom in a magnetic field at scaled energy E = —0.1. s^™ and |A qm |: 
Actions and absolute values of amplitudes obtained by harmonic inversion of the 
quantum spectrum. Sp and \A d \: Actions and absolute values of amplitudes ob- 
tained by periodic orbit theory. 



Re s^/2tt 


Im s^/2tt 


|^qm| 


4/27T 


|A«*| 


0.67746349 


0.00000000 


0.29925998 


0.67746283 


0.29929657 


1.09456040 


0.00003414 


1.07896385 


1.09457049 


1.14785331 


1.11461451 


-0.00000723 


0.51757197 


1.11457036 


0.51991351 


1.35492782 


0.00000089 


0.28203634 


1.35492566 


0.28236677 


1.56500143 


0.00000054 


0.42062528 


1.56499821 


0.42053218 


1.69802779 


-0.00000017 


0.39860976 


1.69802585 


0.39857924 


1.79106306 


0.00000075 


0.34636062 


1.79106067 


0.34672899 


1.87643070 


0.00000051 


0.34352848 


1.87642962 


0.34375997 


1.93352314 


0.00000177 


0.31816675 


1.93352213 


0.31888913 


1.99328138 


0.00000352 


0.32173573 


1.99328294 


0.32238331 


2.03199050 


-0.00017651 


0.31129071 


2.03194819 


0.30888967 


2.07515409 


0.00002170 


0.31251061 


2.07517790 


0.31598033 


2.10239191 


0.00003926 


0.30768786 


2.10234679 


0.30947859 


2.13389568 


0.00007623 


0.33987795 


2.13380984 


0.31999751 


2.15304246 


0.00035849 


0.34510366 


2.15300805 


0.31929721 


2.17505591 


-0.00021399 


0.24582893 


2.17556880 


0.33461616 


2.18921873 


0.00034303 


0.79129103 


2.18914098 


0.84694565 


2.20287712 


-0.00045692 


0.35778573 


2.20439817 


0.36428325 


2.30616383 


-0.00000277 


0.47937694 


2.30615714 


0.47918792 


2.40955313 


0.00000007 


0.35374810 


2.40954609 


0.35398381 


2.58910596 


-0.00002565 


0.57871391 


2.58913149 


0.57621849 


2.60055850 


-0.00002093 


0.94842381 


2.60051951 


0.94785988 


2.62069577 


-0.00001446 


0.39565368 


2.62066666 


0.39496450 


2.70985681 


-0.00000436 


0.09049970 


2.70985132 


0.09040051 


2.76685370 


-0.00008167 


0.70094383 


ghost orbit 


2.84036264 


-0.01311947 


0.00065741 


ghost orbit 


2.87476026 


-0.00973465 


0.00390219 


ghost orbit 


2.95857109 


-0.00000045 


0.31520915 


2.95857293 


0.31514488 
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Table 2 

Hydrogen atom in a magnetic field at scaled energy E = —0.1. Relative deviations 
between the quantum mechanical and classical actions and amplitudes of the three 
recurrence peaks around s/2ir ss 2.6. 



s d /2ir \s^ m - s cl \/s cl |^ qm -^ cl |/|^4 cl | 



2.589131 
2.600520 
2.620667 



9.86 x 10 
1.50 x 10 
1.11 x 10 



0.00433 
0.00060 
0.00174 



Table 3 

Periodic orbit quantities for the three disk scattering system with R = 1, d = 6. 
P: symbolic code; Lp and A^ 1 ' /A^ ': Periods and ratio of zeroth and first order 
amplitudes obtained by harmonic inversion of spectra (hi) and by the classical 
calculations (cl) of Alonso and Gaspard in Ref. [23]. 



p 


rhi 


rcl 
Lp 


\ A (1)/ A (0)\M 


|^(l)/^(0)|cl 





4.000498 


4.000000 


0.3135 


0.31250 


1 


4.268069 


4.267949 


0.5554 


0.56216 


01 


8.316597 


8.316529 


1.0172 


1.01990 


001 


12.321782 


12.321747 


1.3465 


1.35493 


011 


12.580837 


12.580808 


1.5558 


1.55617 


0001 


16.322308 


16.322276 


1.6584 


1.67009 


0011 


16.585261 


16.585243 


1.8514 




0111 


16.849133 


16.849072 


2.1172 


2.12219 


00001 


20.322343 


20.322330 


1.9961 




00011 


20.585725 


20.585690 


2.1682 




00101 


20.638284 


20.638238 


2.4061 




00111 


20.853593 


20.853572 


2.4387 


2.45127 


01011 


20.897413 


20.897369 


2.6127 




01111 


21.117009 


21.116994 


2.6764 
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Table 4 

Non-trivial zeros Wk, multiplicities dk, and error estimate e for the Riemann zeta 
function. 



k 


Re Wk 


Im Wk 


Re dk 


Im dk 


£ 


1 


14.13472514 


4.05E-12 


1.00000011 


-5.07E-08 


3.90E-13 


2 


21.02203964 


-2.23E-12 


1.00000014 


1.62E-07 


9.80E-13 


3 


25.01085758 


1.66E-11 


0.99999975 


-2.64E-07 


5.20E-12 


4 


30.42487613 


-6.88E-11 


0.99999981 


-1.65E-07 


1.90E-12 


5 


32.93506159 


7.62E-11 


1.00000020 


5.94E-08 


7.10E-13 


6 


37.58617816 


1.46E-10 


1.00000034 


5.13E-07 


1.00E-12 


7 


40.91871901 


-3.14E-10 


0.99999856 


1.60E-06 


4.90E-11 


8 


43.32707328 


1.67E-11 


1.00000008 


3.29E-07 


1.90E-12 


9 


48.00515088 


4.35E-11 


0.99999975 


-1.35E-07 


1.40E-12 


10 


49.77383248 


7.02E-11 


1.00000254 


-4.59E-07 


1.10E-10 


11 


52.97032148 


1.92E-10 


1.00000122 


7.31E-07 


6.00E-11 


12 


56.44624770 


-1.30E-10 


0.99999993 


4.51E-07 


5.50E-12 


13 


59.34704400 


5.40E-11 


0.99999954 


2.34E-06 


2.30E-10 


14 


60.83177852 


-3.94E-10 


1.00000014 


1.11E-06 


3.00E-11 


15 


65.11254406 


-4.98E-09 


0.99998010 


-8.30E-06 


2.70E-08 


16 


67.07981053 


-2.05E-10 


0.99999892 


-8.04E-07 


5.30E-11 


17 


69.54640171 


2.51E-11 


0.99999951 


9.45E-07 


6.80E-12 


18 


72.06715767 


-5.74E-10 


0.99999974 


8.63E-06 


4.30E-10 


19 


75.70469070 


3.93E-10 


1.00000082 


1.07E-06 


3.20E-11 


20 


77.14484007 


-2.70E-12 


0.99999979 


1.25E-06 


1.30E-11 


21 


79.33737502 


-3.58E-11 


1.00000086 


3.03E-07 


9.20E-12 


22 


82.91038085 


1.56E-10 


0.99999912 


-8.58E-07 


1.60E-11 


23 


84.73549298 


3.34E-10 


0.99999940 


-7.09E-07 


2.50E-11 


24 


87.42527461 


1.20E-09 


0.99999866 


1.39E-06 


6.70E-11 


25 


88.80911121 


-9.42E-10 


1.00000101 


1.49E-06 


4.80E-11 


26 


92.49189927 


-4.11E-09 


0.99999761 


-1.93E-06 


1.50E-10 


27 


94.65134404 


-7.11E-09 


1.00000520 


-1.27E-06 


6.80E-10 


28 


95.87063426 


8.06E-09 


0.99999001 


-1.15E-05 


5.20E-09 


29 


98.83119422 


-1.78E-11 


0.99999936 


5.70E-07 


3.10E-12 


30 


101.31785101 


4.22E-11 


0.99999969 


-4.73E-07 


4.50E-12 
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Table 5 

Non-trivial zeros Wk, multiplicities dk, and error estimate e for the Riemann zeta 
function. 



k 


Re Wk 


Im Wk 


Re dk 


Im dk 


e 


2532 


3063.43508648 


-1.64E-09 


0.99999901 


1.34E-06 


5.50E-11 


2533 


3065.28655558 


1.15E-09 


1.00000107 


-3.63E-07 


2.40E-11 


2534 


3066.32025039 


-1.66E-10 


1.00000231 


1.00E-06 


1.20E-10 


2535 


3067.07132023 


3.68E-09 


1.00000334 


2.73E-07 


2.20E-10 


2536 


3068.01350133 


-1.51E-09 


1.00000291 


1.46E-06 


2.10E-10 


2537 


3068.98426618 


-5.92E-09 


1.00000205 


2.94E-06 


2.60E-10 


2538 


3069.78290477 


-4.40E-09 


1.00000237 


2.51E-06 


2.40E-10 


2539 


3070.54262154 


-7.71E-10 


1.00000169 


9.57E-07 


7.90E-11 


2540 


3072.00099337 


-6.44E-11 


0.99999908 


2.17E-07 


2.00E-11 


2541 


3073.18523777 


9.17E-11 


0.99999942 


-1.07E-06 


3.00E-11 


2542 


3074.52349428 


6.73E-09 


1.00000391 


-6.51E-07 


3.50E-10 


2543 


3075.03387288 


-1.22E-08 


1.00000117 


-5.69E-06 


7.30E-10 


2544 


3075.83347924 


-3.13E-09 


1.00000013 


-3.86E-06 


3.10E-10 


2545 


3077.42747330 


5.76E-10 


1.00000561 


4.69E-06 


1.10E-09 


2546 


3078.28622690 


1.34E-08 


1.00001283 


1.10E-06 


3.80E-09 


2547 


3078.89737915 


1.61E-09 


1.00000487 


-8.04E-06 


2.10E-09 


2548 


3079.87139464 


1.70E-09 


1.00000275 


-2.32E-06 


3.00E-10 


2549 


3080.85638233 


8.67E-10 


1.00000159 


-4.12E-07 


5.90E-11 


2550 


3082.16316375 


-5.88E-10 


1.00000013 


8.44E-07 


1.70E-11 


2551 


3083.36135798 


8.43E-10 


0.99999923 


3.45E-07 


1.50E-11 


2552 


3084.83845150 


2.72E-09 


1.00000057 


-2.86E-06 


1.80E-10 


2553 


3085.37726898 


-1.37E-08 


0.99999576 


-2.88E-06 


5.50E-10 


2554 


3085.96552225 


6.39E-09 


0.99999667 


1.50E-06 


2.80E-10 


2555 


3087.01881535 


3.46E-11 


0.99999845 


-3.63E-07 


5.20E-11 


2556 


3088.08343703 


-3.89E-10 


0.99999931 


-8.44E-07 


2.40E-11 


2557 


3089.22230894 


-3.31E-10 


1.00000017 


-9.21E-07 


1.80E-11 


2558 


3090.28219490 


2.97E-10 


1.00000069 


-7.17E-07 


2.10E-11 


2559 


3091.15446969 


1.10E-09 


1.00000052 


-6.59E-07 


1.50E-11 


2560 


3092.68766704 


2.25E-09 


1.00000033 


1.45E-06 


5.20E-11 


2561 


3093.18544571 


-2.33E-09 


1.00000168 


-1.50E-07 


6.40E-11 



Table 6 
Trivial zeros and pole of the Riemann zeta function. 



Re Wk 


Im Wk 


Re dk 


Im dk 


£ 


0.00000000 


0.50000000 


-1.00000002 


-4.26E-08 


1.80E-14 


-0.00000060 


-2.49999941 


0.99992487 


-3.66E-05 


1.80E-07 


-0.00129915 


-4.49987911 


1.00069939 


-3.25E-03 


4.40E-05 


-0.09761173 


-6.53286064 


1.07141445 


-1.49E-01 


1.70E-03 
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Table 7 

Weight factors J\f p for the symmetry decomposition of the three disk scattering 
system. 



C?, v 


Al 


A 2 


E 


e 


1 


1 


2 


c 3 ,c% 


1 


1 


-1 


<?v 


1 


-1 
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Table 8 

Semiclassical and exact quantum mechanical resonances for the three disk scattering 
problem (A\ subspace) with R = 1, d = 6. Resonances k hl have been obtained by 
harmonic inversion of the periodic orbit sum, resonances k ce (cycle expansion) and 
the exact quantum resonances k qm have been calculated by Wirzba. 



Re k M 


Im k M 


Re k ce 


Im k ce 


Re fc qm 


Im fc qm 


0.7583139 


-0.1228222 


0.7583139 


-0.1228222 


0.6979958 


-0.0750137 


2.2742786 


-0.1330587 


2.2742786 


-0.1330587 


2.2396014 


-0.1187664 


3.7878768 


-0.1541274 


3.7878768 


-0.1541274 


3.7626868 


-0.1475455 


4.1456898 


-0.6585397 


4.1474774 


-0.6604761 


4.1316606 


-0.6170418 


5.2960678 


-0.1867873 


5.2960678 


-0.1867873 


5.2756666 


-0.1832203 


5.6814976 


-0.5713721 


5.6820274 


-0.5715543 


5.6694976 


-0.5534079 


6.7936365 


-0.2299221 


6.7936365 


-0.2299221 


6.7760661 


-0.2275078 


7.2240580 


-0.4954243 


7.2242175 


-0.4954066 


7.2152706 


-0.4856243 


8.2763906 


-0.2770805 


8.2763906 


-0.2770805 


8.2611376 


-0.2749083 


8.7792134 


-0.4302561 


8.7791917 


-0.4302718 


8.7724709 


-0.4241019 


9.7476329 


-0.3208170 


9.7476329 


-0.3208170 


9.7345075 


-0.3188052 


10.3442257 


-0.3781988 


10.3442254 


-0.3781988 


10.3381881 


-0.3737056 


11.2134778 


-0.3599639 


11.2134779 


-0.3599639 


11.2021099 


-0.3582265 


11.9134496 


-0.3357346 


11.9134496 


-0.3357346 


11.9075971 


-0.3322326 


12.6775319 


-0.3961154 


12.6775320 


-0.3961159 


12.6675941 


-0.3946675 


13.4826489 


-0.2969478 


13.4826489 


-0.2969477 


13.4769269 


-0.2941108 


14.1424136 


-0.4300604 


14.1424117 


-0.4300584 


14.1337039 


-0.4288264 


15.0473050 


-0.2578357 


15.0473050 


-0.2578357 


15.0416935 


-0.2555072 


15.6114431 


-0.4603838 


15.6113293 


-0.4604377 


15.6037211 


-0.4592793 


16.6025598 


-0.2188731 


16.6025599 


-0.2188731 


16.5970551 


-0.2170025 


17.0875557 


-0.4826796 


17.0876372 


-0.4827914 


17.0809957 


-0.4815364 


18.1465009 


-0.1842319 


18.1465009 


-0.1842319 


18.1411380 


-0.1827959 


18.5733904 


-0.4913642 


18.5731865 


-0.4914087 


18.5673580 


-0.4899438 


19.6808375 


-0.1575927 


19.6808375 


-0.1575927 


19.6756560 


-0.1565422 


20.0679755 


-0.4814959 


20.0685648 


-0.4842307 


20.0633947 


-0.4825586 


21.2080634 


-0.1403086 


21.2080634 


-0.1403086 


21.2030727 


-0.1395793 


21.5736413 


-0.4643351 


21.5736471 


-0.4643077 


21.5689872 


-0.4625282 


22.7296581 


-0.1323433 


22.7296581 


-0.1323433 


22.7248396 


-0.1318732 


23.0872357 


-0.4363434 


23.0872405 


-0.4363396 


23.0829816 


-0.4345679 


24.2458779 


-0.1331136 


24.2458779 


-0.1331136 


24.2411967 


-0.1328421 


24.6079881 


-0.4038977 


24.6079949 


-0.4038962 


24.6040788 


-0.4022114 


25.7560497 


-0.1415265 


25.7560497 


-0.1415265 


25.7514701 


-0.1413856 


26.1353684 


-0.3694453 


26.1353709 


-0.3694420 


26.1317670 


-0.3679001 


27.2592433 


-0.1556272 


27.2592433 


-0.1556272 


27.2547515 


-0.1555416 


27.6694386 


-0.3353429 


27.6694384 


-0.3353414 


27.6661128 


-0.3339924 


28.7553526 


-0.1727506 


28.7553526 


-0.1727506 


28.7509680 


-0.1726604 


29.2098235 


-0.3041527 


29.2098230 


-0.3041533 


29.2067213 


-0.3030272 
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Table 9 
Same as Table 8 but for resonances with 120 < Re k < 132. 



Re /c hi 


Im /c hi 


Refe ce 


Im k cc 


Re k* m 


Im fc qm 


120.0966075 


-0.1313240 


120.0966093 


-0.1313207 


120.0956703 


-0.1313419 


120.3615795 


-0.4242281 


120.3624892 


-0.4246464 


120.3617516 


-0.4244457 


120.8941439 


-0.5145734 


120.8959933 


-0.5143028 


120.8949118 


-0.5143053 


121.2644970 


-0.4017485 


121.2643338 


-0.4018042 


121.2633401 


-0.4017384 


121.6158735 


-0.1451572 


121.6158734 


-0.1451572 


121.6149357 


-0.1452113 


121.9157839 


-0.3979469 


121.9158422 


-0.3978279 


121.9151103 


-0.3976081 


122.3933221 


-0.5441615 


122.3962246 


-0.5427722 


122.3951716 


-0.5427575 


122.7533968 


-0.3809547 


122.7533176 


-0.3809838 


122.7523307 


-0.3809219 


123.1345930 


-0.1656546 


123.1345929 


-0.1656547 


123.1336504 


-0.1657315 


123.4680089 


-0.3672776 


123.4679489 


-0.3672571 


123.4672237 


-0.3670296 


123.8916697 


-0.5717434 


123.8969484 


-0.5675655 


123.8959339 


-0.5675273 


124.2423970 


-0.3593830 


124.2423480 


-0.3593825 


124.2413639 


-0.3593262 


124.6520739 


-0.1905948 


124.6520738 


-0.1905950 


124.6511239 


-0.1906812 


125.0191421 


-0.3358345 


125.0191284 


-0.3358584 


125.0184111 


-0.3356413 


125.3924617 


-0.5950322 


125.4000681 


-0.5858492 


125.3990978 


-0.5857818 


125.7306095 


-0.3386874 


125.7305801 


-0.3386762 


125.7295973 


-0.3386268 


126.1681278 


-0.2172657 


126.1681277 


-0.2172664 


126.1671733 


-0.2173496 


126.5700003 


-0.3071799 


126.5700090 


-0.3071895 


126.5692997 


-0.3069994 


126.8986333 


-0.6105834 


126.9067485 


-0.5956839 


126.9058213 


-0.5955889 


127.2175968 


-0.3201029 


127.2175825 


-0.3200895 


127.2166018 


-0.3200477 


127.6830865 


-0.2434140 


127.6830866 


-0.2434159 


127.6821330 


-0.2434894 


128.1211609 


-0.2838964 


128.1211675 


-0.2838934 


128.1204656 


-0.2837384 


128.4113722 


-0.6157741 


128.4171380 


-0.5966404 


128.4162452 


-0.5965218 


128.7033407 


-0.3044266 


128.7033374 


-0.3044145 


128.7023607 


-0.3043815 


129.1973295 


-0.2678886 


129.1973310 


-0.2678924 


129.1963819 


-0.2679547 


129.6731970 


-0.2671784 


129.6731961 


-0.2671733 


129.6725028 


-0.2670549 


129.9292721 


-0.6091832 


129.9304559 


-0.5891942 


129.9295889 


-0.5890492 


130.1879622 


-0.2922354 


130.1879664 


-0.2922273 


130.1869953 


-0.2922052 


130.7109808 


-0.2904524 


130.7109850 


-0.2904576 


130.7100422 


-0.2905099 


131.2271782 


-0.2573647 


131.2271733 


-0.2573650 


131.2264935 


-0.2572819 


131.4488921 


-0.5905439 


131.4451877 


-0.5734808 


131.4443391 


-0.5733029 


131.6713958 


-0.2842971 


131.6714040 


-0.2842951 


131.6704389 


-0.2842886 
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Table 10 

Semiclassical and exact quantum mechanical resonances for the three disk scattering 
problem (A\ subspace) with R = 1, d = 2.5. Resonances k have been obtained by 
harmonic inversion of the periodic orbit sum, resonances k ce (cycle expansion) and 
the exact quantum resonances fc qm have been calculated by Wirzba. 



Re /c hi 


Im k hi 


Re k ce 


Im A; ce 


Re fc qm 


Im fc qm 


4.5811788 


-0.0899874 


4.5811768 


-0.0899911 


4.4692836 


-0.0015711 


7.1436576 


-0.8112421 


7.1441170 


-0.8107256 


7.0917132 


-0.7207876 


13.0000861 


-0.6516382 


13.0000510 


-0.6516077 


12.9503201 


-0.6282388 


17.5688087 


-0.6848798 


17.5699350 


-0.6845598 


17.5042296 


-0.6352583 


18.9210923 


-0.7860616 


18.9266454 


-0.7836853 


18.9254527 


-0.7662881 


26.6126422 


-1.8066664 


26.5430240 


-1.8607912 


26.5316680 


-1.8273425 


27.8887512 


-0.5431988 


27.8881874 


-0.5431833 


27.8577850 


-0.5499256 


30.3885091 


-0.1134324 


30.3884554 


-0.1134504 


30.3528872 


-0.1056653 


32.0961916 


-0.6218604 


32.0966982 


-0.6223710 


32.0693735 


-0.6077432 


33.7905462 


-1.9974557 


33.9721550 


-2.0780146 


33.9323425 


-2.0379172 


36.5066867 


-0.3845559 


36.5066384 


-0.3846400 


36.4822805 


-0.3839172 


39.8138808 


-0.3582253 


39.8139242 


-0.3580251 


39.7859698 


-0.3508678 


42.6585345 


-0.3514612 


42.6556455 


-0.3492807 


42.6312420 


-0.3403609 


44.3271565 


-0.3596874 


44.2457223 


-0.4072838 


44.2295828 


-0.4015643 


45.1114372 


-0.4173349 


45.0602989 


-0.3453843 


45.0400958 


-0.3421851 


48.8432026 


-0.5925064 


48.8417782 


-0.5913485 


48.8203073 


-0.5853389 


51.9171049 


-0.6951700 


51.9146021 


-0.6791551 


51.8983106 


-0.6721309 


53.3788306 


-0.1062924 


53.3766455 


-0.1005604 


53.3584390 


-0.1005951 


60.5907039 


-0.8405420 


60.6204084 


-0.8126455 


60.6075221 


-0.8098089 


62.2000118 


-0.2161506 


62.2004010 


-0.2137105 


62.1832903 


-0.2111191 


64.8371533 


-2.1940372 


64.1528609 


-1.6713984 


64.1457569 


-1.6638733 


65.6804780 


-0.2791339 


65.6804671 


-0.2737728 


65.6638670 


-0.2725830 


67.8612092 


-0.2865945 


67.8688350 


-0.2881904 


67.8515132 


-0.2865564 


69.3238832 


-0.2988268 


69.3443564 


-0.3123706 


69.3334558 


-0.3092481 


71.1163719 


-0.5635882 


71.0822786 


-0.5381921 


71.0672666 


-0.5353411 


74.8566924 


-0.3058255 


74.8552672 


-0.3022484 


74.8405304 


-0.2994128 


77.3118119 


-0.3223758 


77.3193218 


-0.3129301 


77.3088118 


-0.3107104 


78.5807841 


-0.9179862 


78.9236428 


-0.9433778 


78.9042646 


-0.9416087 


80.3395505 


-0.3052621 


80.4173794 


-0.3670163 


80.4002169 


-0.3628903 


83.9450649 


-0.4910890 


83.8740934 


-0.5035157 


83.8631132 


-0.5005356 


85.8041790 


-0.3988778 


85.8000989 


-0.4147646 


85.7918872 


-0.4152898 


88.5979364 


-0.6089526 


88.4703014 


-0.6739426 


88.4561398 


-0.6778186 


93.2377427 


-0.1707857 


93.0234303 


-0.1231810 


93.0113364 


-0.1220480 


94.4188748 


-0.5644403 


94.4284207 


-0.5013008 


94.4351896 


-0.5267138 


97.5331975 


-0.4015302 


97.5576754 


-0.4167467 


97.5478240 


-0.4181066 
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Table 11 
Nearly degenerate semiclassical states of the circle billiard, k : Results from 
EBK-quantization. k : Eigenvalues obtained by harmonic inversion of cross- 
correlated periodic orbit sums. States are labeled by the radial and angular mo- 
mentum quantum numbers (n,m). 



n 


m 


k EBK 


k m 


1 



4 

7 


11.048664 
11.049268 


11.048569 
11.049239 


3 




1 
9 


13.314197 
13.315852 


13.314216 
13.315869 


3 

1 


2 

7 


14.787105 
14.805435 


14.787036 
14.805345 


1 
5 


11 
1 


19.599795 
19.609451 


19.599863 
19.608981 


1 
6 


15 
2 


24.252501 

24.264873 


24.252721 
24.264887 
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Table 12 

The 40 lowest eigenstates of the circle billiard with radius R = 1. n,m: Radial 
and angular momentum quantum numbers; A; EBK : Results from EBK-quantization; 
k(°> and k^ 1 ' : Zeroth and first order semiclassical eigenvalues obtained by harmonic 
inversion of the periodic orbit signal; k cx : Exact eigenvalues, i.e., zeros of the Bessel 
functions J m (kR) = 0. 



n 


m 


k EBK 


k (0) 


k& 


k ex 








2.356194 


2.356230 


2.409288 


2.404826 





1 


3.794440 


3.794440 


3.834267 


3.831706 





2 


5.100386 


5.100382 


5.138118 


5.135622 


1 





5.497787 


5.497816 


5.520550 


5.520078 





3 


6.345186 


6.345182 


6.382709 


6.380162 


1 


1 


6.997002 


6.997006 


7.015881 


7.015587 





4 


7.553060 


7.553055 


7.590990 


7.588342 


1 


2 


8.400144 


8.400145 


8.417503 


8.417244 


2 





8.639380 


8.639394 


8.653878 


8.653728 





5 


8.735670 


8.735672 


8.774213 


8.771484 


1 


3 


9.744628 


9.744627 


9.761274 


9.761023 





6 


9.899671 


9.899660 


9.938954 


9.936110 


2 


1 


10.160928 


10.160949 


10.173568 


10.173468 


1 


4 


11.048664 


11.048635 


11.063791 


11.064709 





7 


11.049268 


11.049228 


11.087943 


11.086370 


2 


2 


11.608251 


11.608254 


11.619919 


11.619841 


3 





11.780972 


11.780993 


11.791599 


11.791534 





8 


12.187316 


12.187302 


12.228037 


12.225092 


1 


5 


12.322723 


12.322721 


12.338847 


12.338604 


2 


3 


13.004166 


13.004167 


13.015272 


13.015201 


3 


1 


13.314197 


13.314192 


13.323418 


13.323692 





9 


13.315852 


13.315782 


13.356645 


13.354300 


1 


6 


13.573465 


13.573464 


13.589544 


13.589290 


2 


4 


14.361846 


14.361846 


14.372606 


14.372537 





10 


14.436391 


14.436375 


14.478531 


14.475501 


3 


2 


14.787105 


14.787091 


14.795970 


14.795952 


1 


7 


14.805435 


14.805457 


14.821595 


14.821269 


4 





14.922565 


14.922572 


14.930938 


14.930918 





11 


15.550089 


15.550084 


15.593060 


15.589848 


2 


5 


15.689703 


15.689701 


15.700239 


15.700174 


1 


8 


16.021889 


16.021888 


16.038034 


16.037774 


3 


3 


16.215041 


16.215047 


16.223499 


16.223466 


4 


1 


16.462981 


16.462982 


16.470648 


16.470630 





12 


16.657857 


16.657846 


16.701442 


16.698250 


2 


6 


16.993489 


16.993486 


17.003884 


17.003820 


1 


9 


17.225257 


17.225252 


17.241482 


17.241220 


3 


4 


17.607830 


17.607831 


17.615994 


17.615966 





13 


17.760424 


17.760386 


17.804708 


17.801435 


4 


2 


17.952638 


17.952662 


17.959859 


17.959819 


5 





18.064158 


18.064201 


18.071125 


18.071064 
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